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FOREWORD 


This report (USCAE 151) is based on a study prepared for the article 
"Perspectives on Hypersonic Viscous Flow Research " in the Annual 
Reviews of Fluid Mechanics. Owing to space limitation, the following 
material in §§6-8 was not included in the article that will appear. The 
original review study is documented here as a whole under a slightly 
revised title, reflecting an expanded scope to include aspects of 
nonequilibrium aerothermodynamics, combustion, and rarefied gas 
dynamics. 



1. INTRODUCTION 


The National Aerospace Plane (NASP) and several other space programs initiated during the 
past decade in the U.S. and abroad (see Williams 1 986, Parks & Waldman 1990, Parkinson & Conchie 
1990, Koelle 1990, Ito et al. 1990, Lozino- Lozinsky & Neiland 1989) have rekindled considerable 
interest in hypersonics. Almost one quarter of a century separates the present from the dynamical era 
of hypersonic flow research in the mid-1950s and early 1 960s, during which critical flow physics 
problems posed by atmospheric reentry were identified and solved while many aspects of aerodynamic 
and aerothermodynamic theories were established. What, then, are the issues and advances in this 
field as perceived in the modern setting? The immense impact of the computer revolution on the 
design concept and analysis strategy, the experience with theJSpace-Transportation-System (Space 
Shuttle) program, as well as advances in material and propulsion technologies since the 1970s should 
all have made the modern research environ ment and pr ogress vastl y different fro m th ose of the 
Sputnik-Apollo era. This article examines'issues and advances in current hypersonic flow research 
f perceived to be of interest in theoretical fluid/gas dynamics. The scope and depth of the review are 
necessarily limited, as is the list of cited references, although the latter turns out to be quite extensive 
owing to the diverse nature of the field. Helpful are two recent texts by Anderson (1989) and Park 
( 1 990 ) which. provide useful background material for the discussion of current issues. [See the 
reviews by Cheng (1990) and Treanor (1991).] 

The nature of this diverse field may perhaps be appreciated by considering simplistically the 
flight Mach number and the Reynolds number Re,*, (or the Knudsen number Kn = M^/Re,*,) as 
two driving parameters which control the high-temperature real-gas properties and the molecular- 
transport processes. A lowering of Re,,,, (increasing Kn) as the vehicle ascends to the more rarefied 
atmosphere brings about nonequilibrium in the internal molecular excitations and flow chemistry, and 
in the translational motion of the particles as well. As with viscous and diffusive processes, they are 
controlled mainly by particle-collision events. The speed and altitude ranges of the Space Shuttle and 
the NASP ascent/descent corridors encompass most such nonequilibrium domains. Thus, apart from 
the fluid dynamic aspect of hypersonic viscous flows, one must address issues of nonequilibrium gas 
dynamics affecting the flows of interest, hence the use of "perspectives" in the article’s title. 

In the present framework, the study of viscous hypersonic flow will face transition problems 
of two kinds which represent, in fact, the two major areas of current research: the turbulence 
transition at the high Re range and, at the other end, the transition to the free-molecule limit. Work 
on f ully developed turbulent boundary/shear layers are outside the scope of this review; some recent 
work applicable to turbulence transition in hypersonic boundary layers will nevertheless be noted. 
Readers may find helpful insight on turbulence modelling and CFD for aerodynamic flows offered 
in recent articles by Chapman (1992), Moin, (1992), Cheng (1989) and Mehta (1990). Towards the 
rarefied-gas regime, there are quite a few Direct Simulation Monte Carlo (D$MC) calculations of 
varying themes to be studied and several issues on con tinuum extension are in need of clarification. 
Works on gaseous radiation and scramjet combustion will be cited only in relation to nonequilibrium 
gas dynamics and CFD studies. A recent article by Tirsky (1993) on hypersonic flow research is 
called to author’s attention. The work presents a perspective quite different from this review and may 
otherwise be considered complementary to the material discussed in §§6 - 8 below. 

2. HYPERSONIC AIRCRAFT AS WAVERIDER 

On hypersonic vehicle design and research, Townend (1991) lists three recurrent themes: 1. 
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replacement of expendable ballistic space launchers with reusable aerospace planes, 2. hypersonic 
airlines, and 3. transatmospheric orbital transfer vehicles. Central to all three is research which aims 
at integrating air-breathing propulsion into an aerodynamic design called the "wa vender". This term 
refers to a concept evolved from Nonweiler’s (1963,1990) study which utilized the streamlines behind 
a known shock wave for generating examples of three-dimensional (3-D) lifting bodies in a supersonic 
flow- -a particular example of which is the caret wing generated from streamlines behind a plane 
shock (see Kuchemann 1978, pp. 74-79, 450-514; Stollery 1990). 

2.1 Waverider as a Generic Design; The Breguet Range 

A great number of recent studies and overviews on waveriders were presented in the 
proceedings of an international waverider symposium (Anderson et al. 1990), where substantial 
improvement in lift-to-drag ratio (L/D) and other aerodynamic features over standard configurations 
are reported. The article by Eggers et al. (1990) extensively reviews the aerodynamic design 
development related to the waverider concept and is itself a valuable document in hypersonic 
aerodynamics. At this juncture, it will be refreshing to recall a discussion by Kuchemann on aircraft 
cruising range and his vision of hypersonic flight. 

Kuchemann (1978, pp. 7-9) anticipated a trend of increasing propulsive efficiency r) p with 
flight speed, and a corresponding decreasing trend in the lift-to-drag ratio of waverider aircraft, so 
that the product J 7 „L/D remains roughly constant — being close to the value x. With this, and the 
provision that the fuel carried is not too small a fraction of the all-up weight, Kuchemann concluded 
from the Breguet range formula that a nonstop flight to the farthest point on the globe is feasible even 
if hydrocarbon fuel is used, irrespective of flight speed. For a Mach-8 Orient Express or the NASP 
X-30 at a comparable speed, the cruise would take about two hours. This conclusion is made explicit 
in Kuchemann's (1978, p.551) "spectrum of aircraft" reproduced here in Fig. 1. It shows the 
maximum ranges of four types of aircrafts designed for cruising at four very different Mach 
numbers, each allowing a two-hour flight time. 

2.2 Waverider Wing Studies 

the viscous corrections Skin friction must be included in the performance analysis of a 
waverider wing. Results of optimization which takes into account skin friction have been referred 
to as "viscous optimized" (Bowcutt et al. 1987, Corda et al. 1988). The viscous optimized waveriders 
obtained are seen to differ considerably in shape depending on (turbulence) transition locations 
assumed in the calculations, signifying the critical need of a reliable transition prediction method (§5). 


BOUNDARY-LAYER DISPLACEMENT, PLOW CHEMISTRY AND HIGH- ALTITUDE EFFECTS Owing to the 

boundary-layer displacement effect at sufficiently low Reynolds numbers, wing loading, skin friction 
and surface heating rate may increase significantly for a thin wing at low incidence. Anderson et al. 
(1992) studied examples of waveriders with a 60-meter chord optimized for this viscous-interaction 
(displacement) effect. The displacement effects on lift and drag are affected little by the optimization 
performed in the study, which nevertheless, alter the waverider planform and its thickness 
distribution drastically. This observation signifies a configuration insensitivity at a given Re and M^, 
which could be translated to a greater degree of freedom for the designers. The effect of air 
chemistry on waverider aerodynamics has also been studied (Anderson et al. 1992) but found to be 
small for the examples considered. 

For higher altitude applications, Anderson et al. (1991) studied waverider wings 5-meter in 
length at altitudes of 80-120 km, corresponding to a Knudsen-number range from 10' 3 to a unit 
order, even though the wing was generated by an inviscid procedure. In their study. Potter’s (1988) 
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Fig* 1 The Spectrum of Aircraft (from Kiichemann 
1978): Range covered by a two-hour flight a s a 
function of cruise Mach number for four different 
aerodynamic designs. The value R/2Rg = 0.5 rep- 
resents the distance to the farthest end of a great 
circle route. 



Fig. 2 Schematic of a hypersonic ramjet (from Billig 
1992). 
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bridging function for empirical correlation of rarefied hypersonic-flow data was considered (also see 
Warr 1979, Wilhite et al. 1985). This and other questions may better be discussed in the context of 
low-density hypersonic flows in §8. 

issues with A sharp leading EDGE A question on the practicality of the waverider design arises, 
which concerns the sharp leading edge inherent to the inviscid solution procedure used. An important 
issue was raised by Nonweiler, namely, whether a genuinely sharp leading edge made of available 
materials can survive the heat flux from hypersonic flight without the aid of active cooling. An 
answer was offered by the "conducting plate" theory and experiment (Nonweiler et al. 1971, 
Nonweiler 1990) which show that solid- body conductivity and radiative cooling can together be 
effective in limiting the temperature on a sharp-edged wing. For a 14° wedge-shaped leading edge 
built from material with conductivity comparable to graphite, the maximum temperature on a 75° 
swept wing at a speed of 6.5 km/sec is not expected to exceed 2000°K, according to the study. 
Recent progress in material research (e.g. Sanzero 1990) could make this passive-cooling approach 
more attractive. A thin/slender configuration with or without a sharp leading edge is apparently 
preferred over a nonslender/blunt shape in the quest for a high L/D. This may be essential for the 
cruise economy as well as cross range capability in a transatmospheric operation (Walberg 1985). 

2.3 Integrated Aerodynamic Design 

The merit of a waverider or any aerodynamic design cannot be assessed without considering 
the constraints placed by the power plant installation, propulsion concept and other details in an 
integrated design, (Kiichemann 1978, Townend 1991). Figure 2 indicates the various parts of the 
external and internal flows of a scram jet engine of a generic design and the need for an integrated 
analysis (Billig 1992). The surface pressure on the ramp would add substantially to the total lift L and 
drag D; the rear portion where the burned gas exits takes the form of a "half nozzle" where the thrust 
T is principally derived, and the pressure also contributes to the lift and pitching moment. As an 
integrated system, one may speak of the net thrust (T-D) available for acceleration. As the scramjet 
vehicle ascends to higher altitudes, the ability to accelerate further depends on the precarious balance 
between the diminishing T and D. 

3. VISCOUS INTERACTION: COMPUTATIONAL METHODS 

The fluid dynamics of hypersonic flows is complicated by the interaction of the boundary 
layer and shear layer with shock waves, leading to flow separation and instability not amenable to 
straightforward analyses. The need for numerical solutions to the Navier-Stokes (NS) or other full 
equation systems has been made apparent in Fig. 2, where significant interaction of boundary layers 
with shock/expansion waves occurs in most regions. Note that in the straight precombustion passage 
(called "isolator") in Fig. 2, a shock train (not shown) must form through wave reflection and viscous 
interaction. Complicated shock-shock interaction patterns can create a supersonic jet impinging on 
the cowl lip, causing an unexpectedly high local heating rate as was first investigated by Edney (1968) 
and later by Holden et al. (1988) and Glass et al. (1989); this has yet to be explained by a viscous 
interaction analyses (cf. Fig. 3 reproduced from Weiting 1990). As a prelude to the discussions of 
the following sections, several major approaches to viscous-flow calculations underlying much of the 
current hypersonic flow studies will be noted; their extension to nonequilibrium flow calculations will 
be discussed in specific applications later in §6. Some of the basic computational procedures for 
compressible viscous flow calculations have been elucidated in texts and monographs (e.g., Anderson 
et al. 1984, Hoffman 1989). The CFD approaches of interest here will be discussed in three categories 
according to the level of approximation for the governing equations. 

3.1 Interacting Boundary Layer Equations 
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Fig. 3 Interaction of incident shock and how shock near 
engine cowl leading edge: Edney’s Type IV super- 
sonic jet interference pattern (sketch reproduced 
from Weiting 1990). 
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Fig. 4 Surface oil flow pattern and pressure contours in 
the symmetry plane and downstream plane com- 
puted by a thin-layer NS code (CFL3D) for a ramp 
on flat plate at M\ ~ 14, Re ~ 10* (Rudy et al 
1991). The analysis resolved discrepancies between 
experiment and earlier 2-D calculations. 


Fig. 5 Normalized over pressure in a compressive free- 
interaction for 

a = (T./T„)"+i = 1, 2, 10, 50, oo (Brown et 
al 1990). Note the plateau pressure depends on 
surface temperature, Mach and Reynolds numbers 
through <r and T, defined in the text. 
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Unlike in standard boundary layer theory, the streamwise pressure gradients in the boundary- 
layer equations are not given a priori and are determined in a system coupled to the Euler equations 
governing the outer flows. The steady-state version of the equations allows upstream influence and 
flow separation and is comparable to a composite form of the PDEs in triple-deck theory (e.g., Werle 
& Vasta 1974; Burggraf et al. 1979; also Anderson et al. 1984). 

3.2 Parabolized Navier-Stokes ( PNS ) 

Even if the viscous interaction is strong, the thin viscous layer permits the deletion of most 
streamwise partial derivative terms in the viscous, heat-conducting, and diffusive parts of the NS 
equations. For steady-state applications, the reduced nonlinear equations give the appearance of 
PDEs of parabolic type, hence the name PNS. With a hypersonic free stream, the upstream influence 
can be omitted in the less viscous part of the external flow, and the PNS system may then be 
integrated simply as an initial boundary-value problem by marching in the downstream direction. 
Documented results show excellent capability of the methods in shock-capturing, and in describing 
the interactions on the global scale as well as capturing cross-flow separations. As methods for 
analyzing viscous interaction, however, a short coming of these marching procedures is their 
preclusion of the upstream influence of the downstream condition and reverse flows owing to the 
special treatments of the streamwise pressure term needed to suppress the departure solution (Vigneron 
et al. 1978; Anderson et al. 1984, pp. 433-440; for recent PNS implementation in hypersonic flow, see 
Butta et al. 1990; Tannehill et al. 1990; Krawczyk et al. 1989). 


In passing, one notes that the equations in viscous shock-layer theories (Davis 1970, Moss 
1976, Cheng 1963; also Gupta et al. 1992) may also be regarded as the simpler versions of the PNS 
using shock fitting. They may best be discussed in the context of the continuum extension to the 
rarefied gas dynamic regime in §8, inasmuch as a considerable number of their applications have been 
made during the last decade in studies comparing the continuum and particle-simulation models (Moss 
et al. 1987, Moss & Bird 1985, Cheng et al. 1989,1990,1991). 

3.3. Iterative /Time- Accurate Navier-Stokes 

thin-layer navier-stokes To retain the upstream -influence capability, one may restore the time- 
dependent terms to the PNS equations, solving them as a time marching, initial boundary-value 
problem, or apply an iterative procedure to the (steady) PNS equations with suitable outflow boundary 
conditions. This version is often referred to as the thin-layer NS and is believed to be basic to several 
codes in current use: ARC3D (Pulliam & Steger 1980), F3D (Ying et al. 1986), NS3D (Blottner 1990), 
and also CFL3D (Vasta et al. 1989). A space-marching iterative procedure called the "supra- 
characteristic method" (Stookesberry & Tannehill 1986) also belongs to this class. In the first two 
codes, flux-vector splitting with upwind differencing (Steger & Warming 1980; van Leer 1982) or 
similar techniques are used, these prove to be robust in many shock-capturing calculations (Roe 1 986. 
van Leer et al. 1987). The F3D code has been further developed and used successfully in 3-D 
hypersonic flow analyses (Ryan et al. 1990), nonequilibrium hydrogen-air reaction (Lee & Deiwert 
1990) and the hypersonic flow through an expansion slot in a 3-D ramp (Hung & Barth 1990). TVD 
(Total Variation Diminishing) and similar schemes were used in many of these works to enhance the 
shock-capturing capability (cf. Yee 1987). Compared to the full NS calculations, the thin-layer 
version may represent a considerable saving in computer resources and programming effort; global 
convergence to the steady state can be accelerated with some zonal strategy, as seen from the cited 
examples. 

full navier-stokes calculations Similar remarks apply to the full NS equations and calculations. 
Among the full 3-D NS codes used in current hypersonic flow analyses is the LAURA (Langley 
Aerothermodynamic Upwind Relaxation Algorithm) designed for finite-volume formulation (see 



5 


Gnoffo 1989); LAURA adapts Roe’s averaging to flux components across cell boundaries for the 
convective terms (see Roe 1986) and Harten’s (1983) symmetric TVD scheme (see also Yee 1987). 
Unique in this relaxation procedure is the point-implicit strategy, which is believed to render the 
procedure stable for an arbitrary Courant number without the need of solving large, block-tridiagonal 
matrix equations. The LAURA code has been applied to a variety of nonequilibrium 
aerothermodynamic problems as well as rarefied, hypersonic flow studies (cf. Gnoffo 1990, 
Greendyke et al. 1992), to which we shall return later for comments. 

The procedure of MacCormack’s explicit, time-split, predictor- corrector method 
(MacCormack & Baldwin 1975, Hung & MacCormack 1975) and the implicit version (MacCormack 
1982) solve equations in finite- volume conservation-law form and are supposedly second-order 
accurate in space and time, as elucidated in the text of Anderson et al. (1984). Application of the 
explicit 2-D version by Hung & MacCormack (1975) to flow past a flat plate with a compression 
corner at M w * 14.1, Re,*, * 1.04x10®, agree quite well with experimentally measured surface pressure, 
heating rate, and skin friction (Holden & Moselle 1969) for ramp angles a * 0-18°. The same set of 
experimental data was also compared with solutions by the supra-characteristics method in 
Stookesberry & Tannehill (1986). However, the comparison for a = 24° was not satisfactory in either 
study for the reason be noted shortly. A 3-D version of a similar procedure was successfully applied 
to a complete reentry configuration at Mach 6 by Shang & Scherr (1986), assuming 7 * 1.40 and a 
Baldwin-Lomax (1978) turbulence model (cf. Anderson 1989, pp. 353-359). In the implicit version, 
a stage is added to each of the predictor and corrector steps, where an approximately factorized time- 
dependent operator is applied to implicitly update the unknowns by simply inverting bidiagonal 
matrices as was done in the explicit version. Whereas the procedure is unconditionally stable for 
unbounded time steps At according to linear model analyses, the products pAt/ e( Ax) 2 and ^At/Q(Ay) 2 
are required to be bounded to maintain accuracy. One observes, however, that the latter requirements 
are similar to, and almost as restrictive as, the stability condition for an explicit method applied to 
a diffusion/heat equation, and that the density g in these products can cause severe problems in 
rarefied hypersonic flow applications. 

The second-order accurate, compressible NS solver recently proposed by MacCormack & 
Candler (1989) is virtually a relaxation procedure and appears to be extremely promising for 3-D 
applications, according to MacCormack (1990). Type-dependent procedures have been effectively 
implemented according to the flux- vector splitting algorithm in solving equations in conservation- law 
form (Steger & Warming 1980). The splitting , applied mainly to the streamwise flux, is believed to 
help relaxation convergence by virtue of the increased weight of the diagonal elements in the block 
tridiagonal matrix of the difference equations. In MacCormack & Candler’s (1989) procedure, the 
Gauss-Seidel line relaxation is adopted to solve the mfactored matrix equations, thereby avoiding 
unwarranted errors from the approximate factorization, which slows down convergence. The new 
procedure allows large time steps, and the calculation can be performed on a common work station. 
It is unclear if convergence acceleration routines (e.g. Cheung etal. 1991) would be helpful in further 
enhancing the method’s performance. Recently, this procedure has been adapted to solve the (full) 
Burnett (1936) equations for rarefied hypersonic flows (Zhong etal. 1991 a,b), for which issues with 
boundary conditions and solution uniqueness remain unresolved (see §§8.5, 8.6 below). 

importance of s-D influence We return now to the comparison of Hung & MacCormack’s (1975) 
calculation with Holden & Mosselle’s (1969) measurement for the case with the ramp angle a = 24° 
mentioned earlier. In this case, Rudy et al. (1991) used the thin-layer CFL 3D (Vasta et al. 1989) code 
to demonstrate that spanwise (global) 3-D effects can resolve all the noticeable discrepancies. The 
computed surface oil flow and pressure contours in the symmetry plane and a downstream plane are 
reproduced in Fig. 4. This provides perhaps an excellent example in which 3-D computation has 
proven crucial in settling a fluid dynamic issue which would have been perceived as being 2-D in 
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origin. Among other computer programs currently being used in viscous hypersonic flow study are 
those of Edwards & Flores (1990), Thomas & Neier (1990), and Liu & Jameson (1992). 

4. VISCOUS INTERACTION: THEORETICAL DEVELOPMENT 

We turn next to the development of viscous interaction theory in this section and later to the 
investigation of the related instability problem of hypersonic boundary layers in §5. 

4.1 Viscous Interaction on the Triple-Deck Scales 

Significant global interaction of a laminar boundary layer with an external hypersonic flow 
(M x »l) has been the subject of extensive investigation in the past (e.g. Hayes & Probstein 1959, 
Moore 1964, Cox & Crabtree 1965, and the review by Mikhailov et al. 1971). There is yet another 
more universal and important interactive feature of a boundary layer occurring on a much shorter 
scale noted earlier by Lighthill (1953), which permits upstream influence and separation and became 
the focus of a vast number of theoretical studies two decades later (see reviews by Stewartson 1974, 
1981; Smith 1982, 1986; Sychev 1987). Central to all the recent work is the triple-deck theory which 
stipulates a three-tier stucture made up of lower, main and upper decks, with the streamwise scale 
short enough that a small self-induced pressure rise is sufficient to provoke flow reversal and 
separation. 

4.2 Triple-Deck Theory Applied to Hypersonic Flow 


The basic parameter controlling the triple-deck structure for a locally supersonic external flow 
can be written for the present purpose as (Stewartson 1974) 



where T is a function of wall temperature and wall shear immediately upstream of the interaction 
zone, and the product inside the large bracket is simply the Lees-Stewartson global- interaction 
parameter 



familiar from the classical theory. In the preceding, the subscript "1" refers to the condition 
immediately upstream of the triple deck and the constant C is the Chapman-Rubesin coefficient 
where the asterisks refer to the reference temperature of the hypersonic boundary layer. 
It is apparent from (4.1) that x remains to be the important parameter controlling the viscous 
interaction on both global and triple-deck scales. Characterizing the theory for this flow structure 
are the orders of magnitudes of the thickness ratios of the lower, the main and the upper decks, and 
also the normalized pressure and streamwise-velocity perturbations, which are representable, 
respectively, as 

e 5 , c 4 , e s , e 2 e (4.3) 

(The streamwise length scale A for the triple deck is the same as that of the upper deck.) This version 
of the theory is to be referred to as the standard version and requires the e in (4.1) to be 
asymptotically small, and is clearly inapplicable to a regime where x is not small. There is however 
a Newtonian version of this approach which considers (7 -l)/2 being asymptotically small in addition 
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to M x being large, while allowing an unbounded x (§4.4). The upstream influence through the lower 
deck may be best seen from the formulation of Rizzetta et al. (1978) using the shear r = 3u/9y as a 
dependent variable, in which a Neumann boundary condition for r at the wall (y = 0), after 
eliminating the pressure gradient, is 

dr _ f°°( ,w 

di-dT'h (r ) v ( 4 - 4 ) 

where the second x-derivative makes the elliptic nature of the problem apparent. 

Among the examples (see Stewartson 1974; Smith 1982,1986) is the free-inter action solution 
which is an eigen/departure solution that is admissible if provoked. The latter leads to separation 
and flow reversal in the lower deck, and reaches a pressure plateau downstream; it represents 
physically the precursor at the head of a large recirculation region. For a ramp angle in a suitably 
small range, solutions with recirculation and reattachment on the ramp downstream were obtained by 
Rizzeta et al. (1978). 

4.3 Is Departure Solution Admissible at Large x? 

The foregoing discussion would suggest that departure solutions of the triple-deck theory are 
unlikely at large \ (strong global interaction). A classical example of global interaction at an 
unbounded x is that of an aligned flat plate, for which the self-similar solution at a uniform wall 
temperature yields a self-induced pressure p/p^ proportional to x (Stewartson 1955, Hayes & 
Probstein 1959). Neiland (1970) found, however, that an indeterminancy exists for an expansion of 
this solution in descending powers of x, i-e-, 

JL^Dooxl l+-+a 1 x- an +-] (4.5) 

Poo 

where for a certain exponent n, the constant a a cannot be determined. The finding suggests an 
upstream influence excluded by the solution procedure. Using a tangent-wedge pressure formula, 
and assuming a unit Prandtl number and an insulated wall, Neiland found n = 50.6. This value was 
confirmed subsequently in the analysis of Werle’s et al. (1973) analysis which considers a wide range 
of wall temperature, and in Brown & Stewartson’s (1975) investigation where the eigen solution was 
found to be insensitive to the approximation made on the outer flow. There, the exponent n was 
shown to be a function of the specific-heat ratio 7 and of the wall-to-stagnation temperature ratio. 
These features may nevertheless be reconciled with the triple-deck formalism discussed below. 

4.4 The Theory for 7-*l 

The impasse in the triple-deck theory posed by large x is overcome by the theory of Brown 
et al. (1975) based on small (7-l)/2 and high Mj 2 , which could be called a Newtonian theory (Hayes 
& Probstein 1959) but which is less restrictive than the latter since the assumption of a strong shock 
is not strictly required. Let e, e p and A gauge the orders of magnitude of the velocity and pressure 
perturbations, and the streamwise length scale of the triple deck, respectively. These can be expressed 
in this case for a nonvanishing x explicitly as 

£ = (7 - 1 )(T w /T 0 ) 2 C X - 1/2 , e - (7-1 )(T W /T 0 ) 6 , A - (7- 1 ) 3 (T W /T C ) 6 , (4.6) 

showing that a triple-deck structure is possible for 7— >1. They also suggest that wall cooling should 
make the theory much more accurate. Using the tangent- wedge approximation, the crucial pressure- 
displacement relation in Brown et al. (1975) needed for closure of the interaction problem can be 
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written as 

PP=^7(- A+ P) (4.7) 

where A is a displacement due to the lower deck, and n is a constant of the order (7- 1 )(T W /T 0 ) 6 x 2 in 
the case of small x • This may be compared with p * -dA/dx in the (standard) supersonic triple deck 
theory. Computational studies with this version of the theory have been made for a compressive ramp 
(Rizzetta et al. 1978) and for free interaction (Gajjar & Smith 1983). 

On the other hand, for finite x and small (7 - l)/2, (4.7) leads to 

+ (4.8) 

which implies that the boundary- layer outer edge, hence the flow in the upper deck, is little affected 
by the interaction. In this connection, one may examine whether the triple-deck result can be 
reconciled with Nieland's algebraic eigensolution (4.5). The latter may now be interpreted as 

ajx' 2 " - ~ bjexp.(n (4.9) 

Note that x a x" 1 ^ 2 and that the triple deck is centered at x r Now the free- interaction solution in 
the theory of Brown’s et al. (1975) gives a pressure precursor of the same form as (4.9) with the 
exponent n being identified as n »■ (0.8273)/A “0[(7- 1 )* S (T W /T 0 )' 6 ], which is indeed a large number 
for the 7 and T w /T 0 of interest, as was anticipated. The Newtonian version of the analysis (Brown 
et al. 1975) remains to be completed with the inclusion of the centrifugal correction in the Busemann 
pressure formula; this is expected to alter substantially the pressure-displacement relation (4.7). 


4.5 Critical Influence of Wall Cooling 

For hypersonic flight applications, theory and analysis must take into consideration the effect 
of a low wall-to-stagnation temperature (T w /T 0 «l). It may be noted that the assumption T w /T 0 » 
0(1) is implicit in the standard theory, and the wall temperature need not fall too far below the 
stagnation/recovery level before a significant departure from the standard theory can occur, as the 
following will confirm. 


The analysis of Brown’s et al. (1990) on the triple deck for small x identifies a critical wall 
temperature level T w : , , 

— r - ** 1 4«A3 




k-'«( 2 

yj 

\T- 1 

) 


(4.10) 


where A is a normalized undisturbed wall shear (equal to 0.332 for an aligned flat plate), and w is the 
exponent in the viscosity-temperature relation /i a T"; the Newtonian factor 2/(7- 1) is included to 
indicate its influence but the limit 7 -» 1 was Ml taken. Depending on the ratio T w /T w \ three 
distinct wall-temperature ranges exist: 


(i) Supercritical: T W »T W *’ (ii) Transcritical: T w = 0(T W *), (iii) Subcritical: T W «T W (4.11) 


For the supercritical and transcritical ranges, the set of scale factors e, e p and A is not basically 
different from that of the standard theory 

3 


= A- 2 ^3^i, A = A- 5 7 (^) 


V 1 ’ 


(4.12) 


where 



(4.13) 
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The relation between the pressure rise and the lower-deck displacement for the ranges (i) and (ii) can 
be reduced to 

p i = —j£r( A + vPi) (4.14) 

where v = kv lt and k is a constant of order unity determined by the boundary- layer profiles just 
upstream of the triple deck, independent of Xi- The term in (4.14) is absent from the standard 
theory, and represents a transcritical (cold-wall) effect. For the subcritical range (s w <s w ), the 
gauging parameters of (4.12) must change in order to remain small, to keep the reduced PDE in 
canonical form and to avoid degeneracy in the P-A relation. This is accomplished simply by 
replacing v 1 therein by 1/k, and (4.14) changes over for the subcritical case to 

±(P 1 + A)=-y r >P l (4.15) 

Interestingly, the relative scales of the triple deck, i.e. e and A, no longer depend on the Reynolds and 
Mach numbers in this case and vanish with s w . 

Figure 5 reproduces the results of Brown et al. (1990) for the overpressure in a free interaction 
for a a i/j 4 in the range of 1 < a < oo, with the origin of £' located at the separation point. The 
existence of the transcritical and subcritical s w -ranges was anticipated in Neiland (1990, private 
communication), the length scales therein differ however from those in Brown et al. (1990). The 
analysis of Brown et al. shows clearly the drastic reduction in the triple-deck length scales hence in 
the extent of the upstream influence as s w vanishes. This means that laminar separation can occur but 
becomes more abrupt under a strong cooling. 

5. BOUNDARY-LAYER INSTABILITY AND TRANSITION STUDIES 

Many investigations of flow instability and turbulence transition in hypersonic boundary 
layers have been undertaken recently. The development is helped substantially, perhaps, by the 
sequence of analyses on compressible boundary-layer instability made decades earlier by Mack and 
others (see reviews by Mack 1984,1987a,b; Roshotko 1976). 

Some caution should be exercised at this juncture on the use of the viscosity-temperature (p- 
T) relation in extending the stability analysis to high- temperature real-gas flow, apart from other 
more obvious considerations. According to a recent study (Kang & Kune 1991), for example, the 
viscosity of dissociating iodine at T ■ 1000-2000°K will have a negative slope in the p-T relation, i.e. 
d/i/dT<0; similar properties may occur in other dissociating/ionizing gases and their impact on the 
stability analysis need to be ascertained. The other aspect in need of caution is the assumption of 
translational equilibrium in certain stability and transition calculations, where the combination of low 
Re and high Mj makes the gas-rarefaction effect important. Take for example, a hypersonic 
boundary layer on a slender/thin body, which may have a boundary-layer thickness S of 2% the global 
scale L, or larger; in this case, it can be shown that the local mean free path is of the of S. 

5.1 Parallel-Flow Instability Applied to Compressible Boundary Layers 

Lees & Lin (1946) extended the viscous (Tollmien-Schlichting waves) and inviscid (Rayleigh 
theorem) results of parallel-flow instability to the compressible case. They noted that the condition 
EKgDU) = 0 (with D s d/dy) signifies a maximum angular momentum and plays the same role in 
compressible theory as does D 2 U = 0 (an inflection point) in incompressible theory. Unlike in the 
incompressible case, this generalized inflection point can be found at some U = U, in the compressible 
boundary layer on a flat plate, and therefore neutrally stable waves with phase velocity c = U g can 
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exist. Lees & Lin limited their consideration to 2-D subsonic relative waves, i.e. |U(y)-c|<a(y). This 
rules out the "supersonic relative waves" with |U(y)-c|>a(y) and the possibility that, in a supersonic 
boundary layer, the TS-type waves are most amplified at some oblique (wave) angles. Allowance of 
supersonic relative waves would render possible acoustic wave propagation and reflection within the 
boundary layer, admitting a sequence of higher modes for each phase velocity, as Mack 
(1984, 1987a, b) subsequently found. [The modes are designated/ordered by a number "n" according 
to the sign changes (zero crossings) occurring in the pressure profile.] The second mode turns out to 
be the most unstable for flat plates and slender cones at high Re (inviscid) and also for all Re 
(viscous) at M!>4, as was confirmed by subsequent experiments (discussed below). The neutral 
stability waves, both inflectional and noninflectional, are significant (as they are in classical theory) 
in that they identify with parts of the boundaries delimiting the instability/stability domains of 
interest. Figures 6,7 (reproduced from Mack 1985) present these curves of neutral instability in the 
domain of wave number a and Mach number M x for 2-D inflectional and noninflectional waves, 
respectively. The calculations were made for an insulated flat plate. Note that a sequence of 
noninflectional waves of neutral stability can exist for each c in the entire range Uj < c < Uj + a x , but 
the results for c = Uj shown in Fig. 6 are more important since each curve therein forms a part of the 
boundary for some genuinely unstable domain. 

Among several peculiar features of Figs. 6,7 are the similarity of the two graphs in trends at 
high and low a, and the drastic slope change together with what appears to be a mode-switching 
behavior in Fig. 6, to be delineated in §5.3 below. One unique feature of a boundary layer with high 
Mj is the progressive movement of the generalized inflection point towards the boundary- layer outer 
edge as Mj increases. Thus at high M 1( this location falls inside the "edge layer" (Bush & Cross 1967, 
Lee & Cheng 1969) and the stability analysis must take into consideration the appropriate p-T law. 
An adverse effect of wall cooling must be noted. At Mj = 10, Mack’s (1985) calculations revealed that 
the temporal amplification rates of the second, third and fourth modes at T w /T 0 ~ 0.05 are almost 
twice the corresponding rates for an insulated wall. This was confirmed experimentally at least for 
the second mode. 

5.2 Experimental Studies of Hypersonic Boundary Layer Transition 

There have been primarily three sets of experimental studies on hypersonic boundary layer 
instabilities at Mj - 4.5-8. 5 reported in Kendall (1975), Demetriades (1978) and Stetson et al. (1988). 
These focused on flat plates and cones in wind tunnels and employed hot-wire anemometer 
techniques. Kendall’s experiments confirmed the existence of the second mode and its dominance 
in a hypersonic boundary layer, Demetriades verified Mack’s findings on the adverse wall cooling 
effect on the second mode, and Stetson et al. investigated tip-bluntness, wall cooling, and other 
effects on slender cones. The latter studies and related works are comprehensively reviewed in 
Stetson & Kimmel (1992) who also noted the existence of a harmonic of the second mode unaccounted 
for by the theory. Figure 8, reproduced in part from Malik et al. (1990) shows good agreement of 
Stetson’s cone data at M w * 8 with Mack’s calculation for the corresponding outer-edge Mach number 
M x = 6.8 in the second-mode frequencies near the maximum growth rate. The noticeable difference 
in the magnitude of the peak growth rate was believed to be caused by the inadequate accuracy of 
the mean flow represented by the boundary-layer solution, but the results based on the PNS generated 
mean flow in the study of Malik et al. (1990) were still far from the mark (cf. dashes and dash-dots 
in Fig. 8). However, a more recent analysis by Simen & Dallmann (1992) produces growth rates 
(reproduced in dots) rather close to the measurements near the peak, attributed principally to the 
merit of a version of thin-layer NS used in the mean-flow analysis. 

5.3 Asymptotic Properties at High 

near-mode crossing Mack’s result shown in Fig. 6 indicates the existence of a segment on the 
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Fig. 6 Wave number as function of outer-edge Mach 
number of the neutral inflectional instability modes 
admissible to the inviscid stability equations of a 
compressible boundary layer on an aligned flat plate 
(from Mack 1984). 



Fig. 7 Wave number as function of outer-edge Mach 
number of the neutral aon-tn/feciiona/ instability 
modes admissible to the inviscid stability equations 
of a compressible boundary layer on an aligned flat 
plate (from Mack 1984). 


Boundary Layar ' 

PNS blunt coo* . (Malik et. at 1990) 

— — PNS iWp 

*— Thm-layer NS (Sunn k Dallmann 1992), 

A Experimental data (Steteoa 1963) 



Fig. 8 Comparison of amplification rates of the second 
mode as function of frequency using four differ- 
ent basic-state equations with experimental data for 
boundary layer on a cone at Mi = 6.8, Rti = 1731 
(Results and data from Malik et al 1990, Simen k 
Dallman 1992, Stetson 1983). 
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neutral curve of each (inflectional) mode where the slope da/dMj becomes positive, while the slopes 
are negative elsewhere and on the noninflectional neutral curves in Fig. 7. This segment has been 
alluded to as the "vorticity mode" (cf. Mack 1984), which will be named hereafter the vorticity 
submode (of each inflectional neutral mode). The parts on the neutral inflectional curves with 
negative da/dMj will be called acoustic submode and those of the noninflectional neutral waves are 
called acoustic modes. Another feature noticeable from comparing Fig. 6 to Fig. 7 is the close 
proximity of the acoustic submode of an inflectional nth mode (in Fig. 6) to the acoustic 
noninflectional (n- 1 )th mode or nth mode of Fig. 7, depending on whether the acoustic submode is 
to the left or to the right of the vorticity submode. Furthermore, the segments of vorticity submodes 
tend to form a continuous curve at high M x , as suggested by Fig. 6. 

This and the feature of near-mode crossing of the acoustic submodes noted above is best 
explained by Smith & Brown’s (1990) asymptotic result, which illustrates the switching from a 
vorticity to an acoustic submode along the nth inflectional neutral curve: 

(“-5 <*•» 

where Tj stands for (In the value of 1.788 was arrived at from 7 - 1.40, and E has a 

magnitude comparable to Mj* 2 exp(-2aMj 2 ). Thus the two distinct submodes are separated by an 
exponentially small amount, and the vanishing of the two factors on the left leads respectively to the 
vorticity and the acoustic submodes. While the vorticity and acoustic submodes differ little in their 
phase velocities at high Mj, i.e. c = U, + OfL^M^ 2 ), they differ substantially in the wave numbers, 
as the foregoing discussion and (5.1) have indicated. Supported by their analysis. Smith & Brown 
(1990) propose to classify instability modes at high M lt including the neutral waves, into two main 
kinds: one is the acoustic mode with wave number and growth rate given, respectively, by 

a = CXnMf 2 ), aCi/Uj = CXM^ 6 iy 1 ), (5.2) 

and the other is the single "vorticity mode" with 

a - CX^ ), aCj/U = CXMf 2 r x ). (5.3) 

The latter is far more unstable and significant than the acoustic mode at high Mach number (also see 
Brown et al. 1991). The result (5.2) was also noted by Cowley & Hall (1990). A vorticity mode 
similar to (5.3) was found in the mixing layer considered by Balsa & Goldstein (1990). 

viscosity-temperature law dependence The foregoing asymptotic study was based on a linear 
viscosity-temperature law. Assuming p a /T, Blackaby et al. (1992) found, instead of (5.2), for the 
acoustic mode: 


a = CKnMj- 3 / 2 ), aCj/Uj - OfM/ 7 / 2 ) (5.4) 

and for the "vorticity mode", instead of (5.3): 

q-CXI), cxc l /\J 1 ■ CXMj' 2 ). (5.5) 

The reasons for the change may be traced partly to the difference in the edge-layer behavior of the 
mean-flow structure for p a T (Lee & Cheng 1969) and for p « T“, w < 1 (Bush & Cross 1967). 

instability at large X Blackaby et al. also made an instability analysis for the flat-plate problem 
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in the strong-interaction regime (x -* <»), assuming p <x /T, and found 

a = (KM/- 3 *), aCj/Uj = CXM^' 1 ) (5.6) 

where A = 9 'y/( 6 ' 7 - 1 ), with which both wave number and amplification rates are seen to increase with 
for all 7 > 1, quite unlike (5.2)-(5.5) which corresponds to x -* 0- This is rather surprising since 
little supporting experimental evidence of instability can be found in the literature in this case. 

the first-mode/ts waves In a triple-deck formalism. Smith (1989) found that the parallel-flow 
assumption at high Mj cannot hold for TS waves even at a rather low Mj - (XRe 1 / 16 ), and that, in 
order to keep the TS waves effectively subsonic for their validity, they must be directed outside the 
wave-Mach-cone V>, i.e., tan -1. Cowley & Hall (1990) studied the influence of a shock on 

the TS waves in the hypersonic boundary layer over a wedge, using a triple-deck approach. To 
simplify the analysis they introduce a special kind of hypersonic Newtonian approximation, in which 
the entire shock-layer thickness becomes narrow enough to be comparable to the upper deck. 
Seddougui et al. (1991) found an adverse wall-cooling influence on the spatial growth rate of the 
(viscous) TS mode, and showed destabilization of the otherwise stable modes. 

nonlinear evolution of the acoustic mode Nonlinear spatial evolution of the unstable waves 
in the acoustic mode was studied by Goldstein & Wundrow (1990). As seen from Eq. (5.2), the 
amplification rate of this mode/submode at high Mj is so weak that even a pressure fluctuation of 
the order Mf^toMj)* 1 suffices to initiate the nonlinear evolution. 

5.4 Hypersonic Boundary-Layer Transition 

Transition prediction requires the identification of the free-stream disturbance field, and 
determination of the boundary/shear-layer response as well as (linear and nonlinear) amplifications 
of these internalized disturbances prior to the breakdown to turbulence. Following Malik et al. 
(1990), the events up to the nonlinear breakdown will be called the stage of "transition onset", to be 
distinguished from the downstream "transitional zone" that follows. 

transition-onset stage The e N method and variants for locating the onset of transition may still 
work at high (though with N being substantially reduced from the magic "9"; see Malik et al. 
1990); their application requires knowledge of the more dominant first and second (if not all) 
instability modes, which provide the (spatial) amplification rate o in 

N = ln(A t /A 0 ) = f ** o(x)dx (5.7) 

where A c is the amplitude of the internalized disturbance at the onset of instability and the subscript 
"t" signifies the end of the transition-onset stage. An alternative to the e N method is to assign an 
amplitude level to A t at the end of the onset stage, instead of assigning a level for N which controls 
the amplitude ratio A t /A 0 . This calls for the computation of A 0 from the free-stream disturbance 
amplitude via (linear) receptivity theory. A more viable method appears to be combining the 
receptivity, linear stability and secondary instability theories [in the sense suggested by Herbert 
(1988)]. Gortler vortices in the boundary layer of a Mach 5 nozzle have been observed by Beckwith 
& Holly (1981); with Gortler vortices as primary disturbances, for example, the secondary instability 
waves can develop amplitudes comparable to the primary level at high Mj (Spall & Malik 1989, Malik 
& Hussaini 1990) as was found at low Mj 

transitional zone modelling In modelling the transition -zone flow, the Reynolds averaged 
Navier-Stokes (RANS) approach deals with an equation system developed from higher-order moments 
of the ewsemWe-averaged NS equations (e.g. Cebeci & Bradshaw 1988) whereas the Large-eddy 
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simulation (LES) models the turbulence in the subgrid scale and deals with numerical solutions to the 
spatially averaged NS equations (e.g. Reynolds 1976, Lesieur 1990). Although the LES is still in a 
developing stage, its prospect as a reliable flow-data source seems high, even for NASP applications 
(Zang et al. 1989). Using meshes and time steps small enough to resolve the Kolmogoroff scale, direct 
numerical simulation (DNS) may need no subgrid modelling and yield data with adequate details for 
some highly idealized, otherwise costly transitional-flow computations. The latter are essential for 
the LES and RANS calibration. There is also an ONERA/CERT version of the RANS which receives 
encouraging support from comparison with DNS results (cf. Malik et al. 1990). Among the RANS 
arsenal in current development is the k-w model (Wilcox 1988,1991) which captures certain 
transitional-zone properties well. 

Distinct from the RANS and LES approaches, and perhaps more appealing, is the application 
of a nonlinear transition theory by Ng et al. (1990) to the analysis of a secondary instability of the 
(primary) second mode on a cylinder at high M 2 in the early stages of the transition (cf. Malik et al. 
1990, Fig. 17 therein). Their results on the Reyolds-stress profile agrees quite well with the DNS data 
and reveals the predominance of a secondary instability in the vicinity of the critical layer attributed 
to a 3-D nonlinear effect. The analysis may explain the "rope-like structure" observed in the vicinity 
of the critical layer in several transition experiments in hypersonic flows (e.g.. Potter & Whitfield 
1969). 

6. NONEQUILIBRIUM AEROTHERMODYNAMICS: MODELLING AND APPLICATIONS 

The works to be examined represent a current research development in high- temperature flow 
physics and may still be considered far from achieving a methology base with unquestioned certainty; 
the research has identified nevertheless several vital issues of modelling at the atomic/molecular levels 
and provided usef ul engineering estimates. 

6.1 Modelling Transition Among Internal States of Atoms <E Molecules 

the master equations A rational way to derive an equation set for the nonequilibrium 
thermodynamics of interest is to seek the time rate of population change of atoms/molecules at a 
specific (energy) state i as the difference between the sum of rates of all collisional (and radiative) 
transitions that populate the state i and the sum of rates that depopulate the state i. Such a system is 
commonly referred to as the Master equations and will furnish the "source term" in a conservation 
equation for N ; (the species population at the state i). The transition/emission rate for each 
collisional/radiative process is to be determined with models/approximations of quantum mechanics 
of the molecular, atomic and electronic interactions, and with the help of the detail balance hypothesis 
(see e.g.. Park 1990, pp. 90-92; Clarke & McChestney 1964, p. 325). 

qss model or electronic STATE Several available computer codes for predicting radiation 
intensities from flows in thermo-chemical nonequilibrium, such as the NEQRAP and NEQAIR 
quoted in Park’s book, were based on a quasi -steady- state (QSS) model of the master equation for the 
(atomic) electronic states, which determines N; in terms of a free-electron temperature T e , and the 
electron-number density N e . The latter is determined separately by a rate equation. In NEQAIR, it 
was assumed that the system could be characterized by three temperatures, T,, T v and T; the 
vibrational temperature T v is shared by all molecules and the translational temperature T is shared 
by all heavy particles. 

6.2 Vibrational Relaxation/ Excitation and Dissociation 

The vibration excitation is considered the main channel of energy transmission to the upper 
level for dissociation and thus controls the dissociation process. For the study of vibrational 
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nonequilibrium, radiative transitions are considered to be much slower than that by collisions and are 
deleted from the master equation governing the vibrational transitions (cf., e.g. Park 1990, p. 97). 

vibration-translation energy exchange Here, the rates of state-to-state transitions are to be 
furnished by modelling collision and excitation processes and may lead to quite different results 
depending on the forms of intermolecular and interatomic potentials and other approximations used. 
If one considers an end-on colinear (aligned, 1-D) collision of a rotationless harmonic-oscillator 
molecule with an atom or a molecule which has a frozen vibrational energy, assuming also an 
exponentially decaying interaction potential as in Landau & Teller’s (1936) theory, one finds that 
transitions can occur only to the neighboring states (v' = v±l) and a deactivating transition rate 
proportional to v. This leads to the familiar relaxation equation for the averaged vibrational-energy 
if, in addition, the vibrational levels are assumed to populate according to a Boltzmann distribution 
at some vibrational temperature T v , namely, 

Af v = r- 1 (r v (T)-£v] (6-1) 

where T y * denotes the average equilibrium vibrational energy at the translational temperature of the 
(heat-bath) particles, T. The temperature and pressure dependence of the vibrational relaxation times 
r for a number of important air-species pairs have been estimated and correlated with experimental 
data at T up to 5000°K by Millikan & White (1963) and others (see Park 1990). 

vibration-vibration energy exchange Schwartz, Slawsky & Herzfeld (1952) considered the end- 
on colinear collision model of a diatomic molecule pair, of which the vibrational energies in both can 
be activated, assuming again no rotational energy exchange and a form of (separable) intermolecular 
potential in this case, as 

U a exp.[-a(r-0 A r A -0 B r B )] ( 6 2 ) 

where r A and r B are the internuclear separations for molecule A and molecule B, respectively, and 
r is the distance between the mass centers of the tw o molecules. According to the SSH analysis 
specialized to a harmonic oscillator model, the contribution to the rate of transition through 
vibrational energy exchange during collisions are found to differ substantially depending on whether 
the total internal energy change AE = (E' A +E' B ) - (E A + E B ) is (nearly) zero or not; the contributions 
from collisions with AE = 0 (resonant case) were seen to be predominant, and thus preferential . (cf. 
Park’s (1990) discussion, p. 61; the resonant transition may not hold for the high vibrational levels.) 

The SSH analysis extended to models with anharmonic oscillators and dissociation yields quite 
different results, made apparent by Sharma’s et a!, (1988) study elucidated in 
Park (1990). There, the nonuniform vibrational-energy spacing were calculated for several 
intermolecular- potential models including the Sorbie-Sorret as well as a 2-term Dunham potentials, 
the multiple-state transitions v ■ v±2 were also considered. We shall return to this version of SSH 
extension shortly. 

In passing, an earlier study on the anharmonicity by Treanor et al. (1968) should be recalled 
[also see Rich & Treanor’s (1970) review]. Limiting the model to the ladder climbing exchange v' = 
v± 1 , and assuming that V-T transitions are generally much slower than V-V transitions, these authors 
noted a QSS-type solution to the V-V dominated Master equation 

N v * N 0 e' irv exp(-E v /kT) ( 6 - 3 > 

where N c , ^ and T may vary slowly with time at rates comparable to the V-T rates. This result can 
be rewritten as 
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N v « N^xpf-vEj/kT^expKvEj-E^/kT] (6 . 4) 

showing at once that for a simple harmonic oscillator the distribution in (6.3) is that of Boltzmann, 
but for an anharmonic oscillator, the last factor in (6.4) furnishes a needed correction. For T<T V 
corresponding to a sudden drop of translational temperature from an equilibrium condition, 7 is 
negative and (6.3) signifies a population inversion. The study indicated a relaxation time much shorter 
than that from the harmonic-oscillator model at the lower-temperature range. Calculations with an 
extended anharmonic version of the SSH theory allow also molecular rotation of 0 2 , NO, CO, OH, 
H 2 and N 2 ; their application to expanding flows have been carried out recently by Park (1992), Ruffin 
& Park (1992), and Roany et al. (1992), confirming the essence of the non-Boltzmann distribution 
(6.4). 

issue on bottleneck in V-V exchange Returning now to Sharma’s et al. (1988) work for the 
collisional V-V exchange model of anharmonic oscillators, the overall transition rate coefficients 
K(v,v+1) and K(v,v+2) as well as K(v,c) were computed up to the vibrational level v = 50 for 
rotationless N 2 at T = 8000°K, T v = 4000°K. Here, the K(v,c) is the rate coefficient for the transition 
from the vibrational state v to the dissociated state unaccounted in Treanor et al. (1968). The case 
considered with T>T V corresponds to a sudden heating of the gas, such as that occurring behind a 
shock. Interestingly, the calculation indicates a vibrational excitation "bottleneck" around v - 20 
where K(v,v+1) as well as a second moment of K(v,v0 have an extremely low minimum (reproduced 
in Park 1990, Figs 2.10, 3.3). Using this set of K(v.v') and K(v,c), the time-dependent master 
equation yields an evolutional solution for N v which is highly non-Boltzmann with three distinct v- 
ranges, and evidently this bottleneck inhibits the transfer of vibrational energy to the upper states. 
With these non-Boltzmann results, Sharma et al. computed the rate of total vibrational population 
removal, i.e., the forward dissociation rate kf, as well as the rate of total vibrational-energy loss, or 
the average removed vibrational energy e v . 

As pointed out by Gonzales & Varghese (1991), these highly interesting results from an extended 
SSH model are affected, however, by uncertainties related to a number of assumptions made in 
Sharma’s et al. (1988) calculation where corrections for several errors and 3-D effects in the original 
SSH analysis (Schwartz & Herzfeld 1954, Tanczos 1956) were not made. An issue was also raised on 
the "distorted wave approximation" implicit in the SSH theory which may not be appropriate for 
transition involving high relative velocities and multi-quantum transitions (Clarke S McChesney 
1964). A more recent study by Landrum & Candler (1991) on vibration-dissociation coupling in N 2 
used a corrected and updated version of the SSH theory, including also contributions from colinear 
collisions of diatoms and atoms, but the important range parameter "a" of the intermolecular potential 
appears to be improperly determined with the Murrell & Sorbie potential which is appropriate only 
for the field between bound atoms, according to Gonzales & Varghese. Thus a demonstration of the 
existence of the bottleneck in question has not been concretely established. It was also pointed out 
that a recent DSMC calculation by Olynick et al. (1990) using inelastic cross-section data 
corresponding to the rates used in Sharman et al, (1988) gives no sign of such a bottleneck. 

alternative semi-classical approaches Transition-rate calculations in this case could be 
improved by using an extension of the semi-classical "N-state method" (Rapp & Kassal 1969) in 
which the dynamics of the colliding pair is modelled classically while the oscillator motion is modelled 
quantum-mechanically. The task is rather computationally intensive, only a few among the many 
needed state-to-state transition rates have been calculated by Gonzales & Varghese (1991, 1992) to- 
date. Another promising semi-classical approach currently being developed by Kune (1991) treats 
the V-T exchange for undissociated diatomic gases based on a multi-center potential model developed 
earlier (Kune 1990). The treatment includes a quantum-mechanical interaction of vibration and 
rotation, retaining however the simplifying features of the 1-D colllisional model. Sample transition 
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probabilities Pjj.j have been computed for j = 1,4 and 20 for N 2 and 0 2 . Comparison with the 
collision numbers for the 1—* 0 transition, Z 10 inferred from exisiting experimental data (after 
adjusting a constant "g" in the theory) show very encouraging agreement in the temperature range T 
- 1000-6000K 0 . Comparison of this approach with experiment and with other theories for transition 
rates at the higher vibrational levels remain to be made. 


diffusion theory The master equations for the vibrational-state population, allowing also 
transitions to the unbounded (dissociated) state, may be approximated by an integral equation. With 
appropriate assumptions, the latter may in turn be reduced to a diffusion-like equation (Keck & 
Carrier 1965) which is amplified here for its theoretical importance in the interpretation of 
computational and experimental results. If the vibrational energy gap E(v+ 1 )-E(v) is small compared 
to kT, the master equation for the diatomic molecules may be expressed in an integral equation form 
[using mostly Park’s (1990) notations] 

N-'i-piV) = C" K(v,v')Wv') - p{\))dy/ + K(v,c)[p A p B - ?(v)] (6.5) 

dt J 0 

where v and v' are the level number in a suitable unit, N x is the number density of the colliding 
partners, p is the number density at v- level normalized by its equilibrium value, i.e. p m N v /N v ,p A 
and p B are the normalized number density of (free) atoms A and B, and v m is the maximum v. Two 
additional requirements are needed in the diffusion theory; they amount to (i) the rate K(v,v') is large 
only in the vicinity of v' = v, (ii) K(v,c) is appreciable only at those upper states near the dissociation 
limit. Note that the typical mid-level vibrational energy is of the order D, and that in most cases of 
interest kT«D. This, together with the assumption required for the integral form of the Master 
equation, means 

AE ■ E(v+l)-E(v)«kT«D (6-6) 

The kernel K(v,v') is sharply peaked under assumption (i), this allows simplifications of (6.5) to 

N t^p( y ) = ^ [m^/<v)| + K(y, c)\p A pB - P(v)] (6.7) 


where M is a transition moment 

Af(v) = r K(v,v + 0^ 
J — oo 


( 6 . 8 ) 


If one excludes those high lying v-states with energy E(v) close to the dissociation limit D, then, by 
virtue of the assumption (ii) appropriate under (6.6), a diffusion equation follows 

= <69) 

The boundary condition at v ■ 0 is obviously dp/d\ ■ 0 since no molecules can cross this boundary. 
The upper boundary condition at v— *v m was furnished by Keck & Carrier as 

M%- = (PAPB-P) f m K(v,c)dv (6.10) 

ov J o 

An ambiguity then appears since the right member of (6.10) is also seen to be the net formation rate 
of molecules from atoms. A more in-depth discussion on the boundary conditions is given in Park’s 
recent (1992) paper. Equation (6.9) would show that a vanishingly small diffusion coefficient M 
occurring at some mid level would inhibit the upward population migration to the dissociation level 
and would signify a "bottleneck" noted earlier, and the M(v) in the example of sudden heating of N 2 
in Sharma’s et al. (1988) study did reveal an extremely low minimum. The latter result would have 
made great theoretical impact if not for the several uncertainties brought out earlier. 
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A rather encouraging application of the diffusion theory has been made recently by Lee (1992) 
who considers a model of energy exchange of vibrational oscillators with an electron heat bath in the 
absence of dissociation. This model has been considered important in the study of highly ionized re- 
entry flow (Lee 1985). Reasonable agreement of the diffusion model with the original discrete master 
equation in the sudden-heating case is demonstrated here; non-Boltzmann distributions were shown 
at an earlier stage of the transient, but a "bottleneck" in K(vjv') was not evident. 

6.3 Applications, Validation and Assessment 

multiple-temperature models and reaction rates The extensive study of the vibrational 
transition models has led invariably to an appreciation of the multiple-temperature concept in 
nonequilibrium aerothermodynamics. These improvements in concept take the forms of allowing 
several independent temperatures, namely, the translational-rotational temperature T shared among 
heavy particles, the vibrational temperature T v shared by all molecules and (depending on the 
requirement) the electron-electronic temperature T e ; these temperatures enter in the empirical rate 
formulas to reflect on the impact of T, T v and/or T e suggested by the model studies. Note must be 
taken to the rate of vibrational energy change (per unit volume) due to collision which may be 
evaluated as a sum contributed by three sources: (i) the Landau-Teller form for V-T and V-V 

transitions modified by a factor depending on T v and T, (ii) vibrational excitation contributed by 
electronic impact (significant for N 2 ), (iii) vibrational energy removal/addition due to 
dissociated/recombination, with the energy per molecule T y based on some preferential dissociation 
model. This sum then enters as a source term in the conservation PDE governing E v . Similarly, the 
rate of the electronic energy change E e owing to the kinetic process is contributed V by seven terms 
resulting from electronic excitation, ionization, ionic recombination, radiation, etc. This rate enters 
again as a source term in the conservation PDE for E t . 

More empiricism is found in the attempt to modify the pre-exponential temperature dependence 
in the forward reaction rate k f , namely, replacing T by the product of an average suggested by Park 

T a = T’T 1 -* (6.11) 

where q varies between 0 and 0.5 in practice. The lack of adequate experimental data to determine 
the constants and exponents of the rates represent an uncertainty in the study which may therefore 
be considered qualitative at best. Another area of ambiguity concerns transport properties of the gas 
mixture and its simplification (Yos 1963, Wilke 1950), which has never been critically tested or 
ascertained at the elevated temperature range of interest. The need for their scrutiny is made 
apparent by the example of iodine considered by Kang & Kune (1990) noted earlier. 

convective and radiative beatings Gupta (1987) studied the thermochemical and radiative 
properties of the sh ock la yer during the Fi re II r eentry s tudy carried out by a Space Shuttle 
experiment (Cauchon 1966); solutions were obtained for NS equations with nonequilibrium chemistry, 
and for viscous shock-layer (VSL) equations but with equilibrium chemistry. Except at altitudes 
higher than 80 km, the Fire II data of radiative intensities, heat fluxes agree reasonably well with the 
equilibrium VSL solutions for a fully catalytic surface as well as the inviscid shock-layer and 
boundary-layer analyses of Sutton (1984). In a code- calibra tion study, Gnoffo (1990) applied the 
LAURA code using park’s 1 1 -species model with T a = /T v T to predict convective heat transfer rate 
of the Fire II test during the early period (corresponding to 1 1.3 km/sec speed and altitude 85-67 km). 
Good agreement with flight test data are found although the heating rates in this range are too small 
compared to those at the lower altitudes to be of significance. 

In a proposed Aeroassisted Flight Experiment (AFE), Hamilton et al. (1991) predicted the 
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stagnation-point heating history of a fully catalytic heat shield with 2.2 meter nose radius using an 
1 1 -species nonequilibrium air chemistry and the VSL approximation. Peak convective heating in this 
case occurs at 78 km altitude at speed 9.2 km/sec and the predicted value reaches 0.5 mega 
watt/(meter) 2 (which is believed to be the limit for the reusable tile on the Space Shuttle). However, 
this VSL analysis ignores the "shock-slip correction" which would result in a 50% reduction of the 
peak value (cf. §8). 

entry into martian atmosphere and return Radiative heating will dominate the surface heat 

flux at considerably higher reentry velocities and at lower altitudes and on a larger body, such as 
during an aerobraking return from Mars. A coupled radiation and ablation injection model of the 
nonequilibrium viscous shock layer was used in a reentry heating study by Gupta et al. (1990) in the 
latter case, assuming a speed 16-8 km/sec, altitude 80-65 km, and a nose radius 3m. The study 
extended an earlier work of Moss (1976) and assessed the impact of using different transport and 
thermodynamic properties, and also different radiation models and, interestingly, showed the 
adequacy of a universal Lewis number 1.40. Over the speed range 12-16 km/sec, the wall heating 
rate was found to vary from 2.5 to 11 MW/m 2 , of which radiation contributes 40-70%; the 
effectiveness of the albative injection of carbon-phenolic is unclear from the study. Entry into the 
Martian atmosphere represents a different aerothermal environment where C0 2 (97%) and N 2 (3%) are 
the main constituents; estimates for entry vehicles with nose radii varying from 1 to 23 meters 
indicate the need of considering a speed range 6-12 km/sec at altitude 30-50 km. Candler (1990) and 
Park et al. (1991) have studied the nonequilibrium nature of this problem. Candler’s NS calculations 
based on an 8-species chemistry without ionization reveals near thermo-chemical equilibrium in most 
parts of the shock layer, attributed to the rather fast C0 2 vibrational relaxation. Applying a 
computationally more efficient VSL analysis to this problem, Gupta et al. (1991) assumed a full 
thermochemical equilibrium but allowed coupling of the shock-layer thermodynamics to the radiative 
cooling. Their results at the lower speed range (6-6.5 km/sec) support Candler’s (1990) observation 
on chemical equilibrium. Convective heat transfers at 8 km/sec speed contributes to 60% of the total 
heat flux for the 1 m.nose radius, and to 23% for the 23 m.nose radius; at 12 km/sec speed, it amounts 
to only 40% for the lm. nose and 2.4% for the 23 m nose. 

electron-number density A series of instrumented probes, called the RAM-C tests, were flown 
at speed 7.65 km/sec and altitudes 7 1 -8 1 km to measure electron-number density around a spere-cone 
with a 0.152 m nose radius (Jones & Cross 1972). Gnoffo (1990) compared his LAURA calculations 
of electron- number density profile using two sets of chemical kinetic rates with the data measured 
by a Langmuir-probe rake. Only qualitative agreement can be achieved (cf. Gnoffo 1990, Fig. 13). 
Whereas, an earlier VSL calculation by Kang et al. (1973) showed a better agreement with the 
measured data. The discrepancy was believed to have resulted from a fully catalytic surface assumed 
in the LAURA calculation inappropriate for a Teflon coated afterbody. 

Candler & MacCormack (1988) assumed a non-catalytic wall in an earlier NS calculation for this 
case and indeed found reasonable agreement in the electron-number density with the RAM-C 
microwave- reflectometer data measured along the cone afterbody (cf. Figs. 2-4 therein). The basic 
numerical procedure of the last work has been noted earlier in §3; two sets of chemical kinetic models 
were tested, one consists of five species N 2 , 0 2 , NO, N and O, the other consists of two more species 
NO + and e\ The electron density in the 5-species set was generated by a special quasi-steady 
approximation which proved to be inadequate. Treatments of nonequilibrium flow models are similar 
to those of Lee (1985) and Park (1990). The program allowed distinct vibrational temperatures for 
different molecular species, which turned out to be very close to one another, thus supporting Park’s 
idea of a common T v shared by all molecules. 

Using a two-temperature version of the LAURA code, Greendyke et al. (1992) carried out a 
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parametric study of the unknown constants/exponents in the nonequilibrium thermochemical models 
and their impact on electron- number density prediction for the AFE experiment. Variations are 
considered in the reference ionization potential of N and selection of rate constants from among Kang 
& Dunn (1972), Park (1987,1990), and their updates, including the option of using Gupta’s et al. 
(1990) equilibrium constant, and also options of imposing limiting cross sections for vibrational 
excitation. A value of the exponent q in (6.1 1) proposed by Hansen (1991), q = 0.1+0.4(T V /T), were 
also considered; it appears to yield only minor changes in this case. The calculations made for an AFE 
model (2.16m nose radius at 78-81 km altitude, velocity 9.7-8.9 km/sec) reveals significant variations 
in the rates of electron-impact ionization with correspondingly large differences in the location and 
magnitude of the peak electron density. The severity of an electron avalanche associated with changes 
in these models was noted. 

more parametric studies Another and perhaps a more extensive parametric study was made 
earlier by Hartung et al. (1991), Mitcheltree (1991) and Hartung (1991) using the two-temperature 
version of th e LAURA code. Hartung et al. s tudied radiat ive emission pr ofiles and radiation spectra 
in the stagnation region for conditions corresponding to a FIRE II flight experiment (al titud e 76-85 
km, speed 11.4 km/sec, nose radius 0.75m). The radiation model used was the Langley Optimized 
Radiative Nonequilibrium (LORAN) code which proves to differ little from Park’s NEQAIR in 
results. The sensitivity study includes again the choice of the exponent q in (6. 1 1 ) for the dissociation 
rate, and of the limiting vibrational relaxation cross section, both of which were shown to be critical. 
Mitcheltree (1991) examines LAURA solutions of translational and vibrational temperatures, 
electron-number density, (^-concentration as well as convective and radiative surface heating rates. 
The flow condition in the study corresponds to an aerobrake of 1-m nose radius at speed 12 km/sec, 
80 km altitude. Rate parameters variations are seen to have an effect, as large as a factor of three, 
on the ionization degree and radiative heating. The results based on Park’s (1987, 1989, 1990) rate 
sets were seen to be affected little by using Gupta’s et al. (1990) equilibrium constant; Hansen’s model 
(6.12) appears to give results virtually identical to that for q = 0 in (6.1 1). 


Hartung (1991) pointed out that the procedure in Park’s NEQAIR code may lead to a negative 
excitation temperature for a bound-free transition which is avoided in LORAN. Comparison of 
predicted emission spectrum in the visible range from LORAN with an AVCO shock tube experiment 
(Allen 1962) at the condition corresponding to the peak radiation point does not however appear to 
be as good as expected, thus requiring further study. On the other hand, there has been better 
agreement of the NEQAIR prediction with the AVCO emission measurement found earlier by Park 
both in equilibrium and nonequilibrium regions (1989; cf. also Park 1990, Figs. 8.24). The comparison 
is reproduced in Figs. 9a,b where the dash curves are results from DSMC models to be discussed in 
§ 8 . 


7. FLUID DYNAMIC PROBLEMS OF HYPERSONIC AIR BREATHING PROPULSION 


7.1 Preliminary Remarks 

Reviews on major developments in the theory and design of hypersonic air-breathing propulsions and 
their key fluid dynamic problems have been cited earlier in §§1-3. The following will discuss some 
recent works on the combustion fluid dynamics. Unlike the flow chemistry of external aerodynamics, 
the overall combustion chemistry of interest is heat releasing (exothermic) and its interaction with the 
fluid dynamics giving rise to many important features of supersonic combustion. Also unlike the flow 
chemistry at low speed, the resident (flow- transit) time in a supersonic combustion may not be long 
compared to the characteristic time of the chemical reaction so that the energy conversion of interest 





Fig. 9 Comparison of calculated and measured spectra for air behind a shock advancing at lOKm/sec. in 
a shock-tube section with initially pressure 0.10 torr (from Park 1989) : (a) near equilibrium region, 
(b) nonequilibri um . Emission spetra prediction based on Bird’s (1987) DSMC-VHS model are shown 
as dashes. 
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is determined by the finite-rate chemistry and may depend also critically on the fluid-mechanic 
mixing process (Ferri 1973, Billig 1992). 

The propulsive performance estimates from the late 1950’s are nearly identical to those that are 
being calculated more than thirty years later, according to Billig (1992); several long recognized 
problems in supersonic combustion research still remain open today. Some of the observations and 
design guidelines offered in Billig’s (1992) review may also prove helpful to newcomers in the field. 
The fuel specific impulse Isp (sec) is related to the propulsive efficiency tj , the vehicle speed V and 
the caloric fuel value H as Isp * Ht? p /V. Studies up to the mid 1970’s indicated that at Mach 6, Isp 
value of H 2 fuel is about 3000 sec, being almost twice the Borane value (1600 sec) for either ramjet 
or scramjet, and, at Mach 10, the scramjet is clearly superior with an estimated Isp slightly more than 
2000 sec for H 2 and 1000 sec for borane (Waltrup et al. 1976, fig. 1). This substantiates a propulsive 
efficiency of rj p of 0.5-0.6 in the Mach-number range 6-10 anticipated in early studies (cf. §2.1). 

In passing, we recall the severe shock interference heating found at the cowl lip of a scramjet 
engine noted in §3 where the high heating rate 30 times the normal stagnation point value has not yet 
been closely predicted to the best of the writer’s knowledge. 

7.2 Scramjet Mixing-Combustion Studies 

Hj - air mixing and combustion Eckland & Northam (1992) studied the effects of certain geometric 
parameters on a combustor performance under conditions corresponding to M^, = 5-7. The basic 
problems analyzed was the mixing and combustion of an H 2 round jet injected (normally) from a hole 
downstream of a step and the effectiveness of the equivalence ratio ER, wall deflection angle 9, and 
the stagnation temperature T 0 in controlling the chemical and mixing efficiencies. (The ER is defined 
as the ratio of injected fuel to that is required for the stoichiometric reaction.) NASA Langley’s 
Reynolds-averaged NS code SPARK was used with finite-rate chemistry of six reacting and one inert 
species. Several measured mixing and chemical efficiencies, defined in the paper, were determined 
and indicate the important roles of chemical kinetics at higher M^. The nozzle hole shape is known 
to significantly alter the jet mixing rate at low speed (Ho & Gutmark 1987), its potential in enhancing 
mixing in supersonic flows remains to be ascertained. 

Another parametric 3-D computational study was made by Kamath et al. (1991) seeking the 
influence of the ER value and flight Mach number M^ on the performance for a conceptual scramjet 
combustor at M m » 3-20, which consists of a ramp and slot fuel injectors. A parabolized version of 
a SHIP3D (NS) code with a k-u> turbulence closure model (Coakley 1983) was used; chemical 
equilibrium was however assumed. Significant increase in the mixing efficiency with ramp injection 
angle, ramp angle and ramp sweep angle were found. 

It is unclear how sensitive is the solution to the turbulence diffusion model used and how an 
alternative 3-D turbulence model may affect the result. Some of these issues on turbulence mixing 
and thrust augmentation techniques applicable to scramjets were discussed in Kumar et al. (1989). 

slot and contour- wall injbctor mixings A number of recent studies focus on fluid dynamic 
means of increasing mixing efficiency in scramjet combustion (and to compensate for the less 
efficient mixing in the case of an axial injector). Of considerable fluid mechanic interest are a series 
of experimental and computational investigations of shock enhanced mixing, exploiting the principle 
of baroclinic (ApxAp) vorticity production (Marble et al. 1987,1990) and its application to contoured 
wall fuel injectors (Waitzel et al. 1991,1992). The small sketches in Fig. 10 describe the basic model 
geometry. To simulate mixing of air with a light gas, helium was injected axially at a speed 
comparable to that of the air flow. The experiments were conducted in a NASA Langley Mach 6 
wind tunnel. The numerical simulations were performed with the SPARK 3D code, assuming laminar 
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Fig. 10 Computed contours of constant helium mass 
fraction downstream of a contoured wall fuel injec- 
tor in a Mach 6 wind tuntoi at Rt ~ 2x 10 7 (Waitzel 
et al 1992). The shear layer at the injector exit was 
assumed in this example to have a boundary-layer 
to injector-height ratio of 0.20. 


Fig. 11 Schlieren photograph of shock deflagration 
around a blunt obstacle at superdetonation speed 
in a stoichiometric H^-air mixture (reproduced from 
Lee Sc Deiwert 1990, photo-record originally from 
Lehr 1972). 




Fig. 12 Variation of thrust coefficient of a ram accel- 
erator with velocity ratio V/Vqj in three different 
propellant mixture (from Hertzberg et al 1991 ) 
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flows. Thus the analysis was expected to capture only large-scale kinematical processes. The key 
player in the baroclinic vorticity production is the oblique recompression shock right behind the 
injector; its intersection with the (non-parallel) high density gradient at the fuel-air interface 
determine the baroclinic growth of the axial vorticity which is expected to enhance mixing. The 
cross-flow patterns from computation (helium mass-fraction contours) and those from wind-tunnel 
measurement (pitot pressure contours at a slightly higher injection velocity) were broadly in accord. 
An example of the computed cross-flow patterns at successive downstream stations are reproduced 
in Fig. 10. The effects of incoming boundary layer height, injector spacing and other parameters 
were investigated. 

The injector performance in the last cited study was observed to be strongly dependent on the 
thickened hypersonic boundary layer on the ramp. Such a pronounced viscous influence was quite 
clearly shown in an earlier computational study of the hypersonic flow through a narrow expansion 
slot by Hung & Barth (1990) who applied a finite volume, TVD version of the thin-layer ARC 3D 
code (cf. §3) to the flow about an injector which was comparable to those sketched in the last figure, 
but only the viscous flow blockage effect in the narrow expansion slot (between neighboring 
injectors) upstream of the mixing zone was investigated. Results obtained for Mach number 5 and 
Re of 10 s indicates a significant (60%) total pressure loss with a surprisingly minor loss in mass flux 
(12%). Wall cooling is shown to reduce the loss and increase the flow expansion substantially. 

7.3 Premixed Shock- Induced Combustion 

Steady and unsteady detonation waves involving a premixed fuel-air combustion behind a shock is 
a classical topic in combustion theory (Karman et al. 1958, Williams 1985). The following will discuss 
problems in recent studies with detonation waves in hypersonic nows. 

PREMIXED COMBUSTION about A shock holder In the validation study for their F3D/Chem 
computer code, Lee 7 Deiwert (1990) analyzed the supersonic combustion of a stoichiometric H 2 -air 
mixture around a blunt obstacle and compared the density field with the photorecord from Lehr’s 
(1972) early experiment. The finite-rate flow chemistry with seven active and one inert species was 
incorporated into an implicit flux -vector splitting code (Ying’s F3D). Applying to an example with 
a 0.75cm nose radius at pressure 0.42 atm and speed 2.6 km/sec (Mj - 6.46), the present method 
reproduced results from a number of existing codes, but cannot reproduce experimental data with any 
of the three selected sets of rate coefficients. Lehr’s Schlieren photograph shown in Fig. 1 1 reveals 
two distinct density demarcations. The inner demarcation was interpreted to be a flame combustion 
front. This front approaches closely the bow shock in the blunter part of the shock delimited by the 
sonic points (right behind the shock), where the reaction is fast enough to provide a relatively thin 
detonation-wave structure. At the sonic point, the Chapman-Jouguet condition will not allow local 
heat release, and beyond the sonic point on the outer part of the shock, heat addition along the 
streamline is again possible but at a much lower rate, owing to the critical temperature dependence 
of the reaction rates. This results in a much longer "ignition delay" which may explain the marked 
departure of this part of the combustion front, which could also be called "deflagration wave". The 
failure to capture this front in computation could be caused by the omission of H0 2 and H 2 0 2 , and 
to ignoring the influence of the nonequilibrium molecular vibration on the reaction rates, familiar 
from the preceding discussion in §6. The authors note specifically that an exceptionally long 
induction time preceding H 2 dissociation may have occurred in the flow. More recently, this 
computational problem is re-examined in Wilson & MacCormack (1992), implementing the 
MacCormack-Candler (1989) procedure discussed earlier with an adaptive grid, which appears to be 
successful in capturing the deflagration front. Only the inviscid model has been considered. 

RAM accelerator and detonation waves The "ram accelerator" is a projectile launcher virtually 
equivalent to a scramjet-in-tube, with propulsive cycles similar to that generating thrust in an 
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airbreathing ram/scramjet (Hertzberg et al. 1988,1991). Fuel premixed with oxygen and inert gas 
components in the tube is heated and ignited by a shock-wave system attached to the high-speed 
projectile which derives thrust from the combustion. The performance of such a device, and the 
thermodynamic states of the combustible propellant gas can be studied by a 1-D analysis as for a 
stationary Chapman-Jouguet detonation wave, but with the inclusion of a thrust term. As a result, 
the (relative) incoming velocity V that allows the maximum heat addition can no longer be the 
Chapman-Jouguet (C-J) velocity, but becomes generally lower for a nonvanishing thrust. Figure 12, 
reproduced from Hertzberg et al. (1991) as a representative sample, compares the 1-D theory (dotted, 
dashed and solid curves) with data inferred from ram-accelerator measurements for three sets of 
premixed propellants (symbols), assuming that the exothermic reactions can carry to their completion. 
The comparison supports the expectation of the 1-D thrust generation analysis, but at speed 
approaching and beyond the C-J value, experimental data show significant thrust generation which 
cannot be explained by the simple theory. Computational efforts were being made to attack the 
problem (Chuck et al. 1991). The thrust generation in the trans-detonation regime and beyond are 
yet to be investigated theoretically. 

computational singular perturbation There have been a few CFD methods which address the 
complex "stiff equations of nonequilibrium combustion and reacting flows using fractional time steps 
and implicit/coupled strategies with a varying degree of success. A recent approach called 
Computational Singular Perturbation (CSP) recognizes the needs to address the nature of the "many- 
time-scale problems" and to design algorithm for deriving computationally time-resolved, simplified 
kinetic models from the otherwise complex reaction system (Lam et acl. 1989, Lam & Goussis 1990). 
The work perceives the concentrations of N species as components of a vector y and their governing 
kinetic equation as 

|y = G(y) (7.1) 

with the global reaction rate G being contributed by the M elementary reactions 

G = J2F i S j (7.2) 

i=i 

where F^ is the reaction rate of the jth elementary reaction and Sj is the stoichiometric (column) 
vector. The key idea of the CSP is to project the M terms in (7.2) into N modes associated with a 
desirable system of N linearly independent basis vectors a s ’s. Each of a s contains a suitably chosen 
time scale, so that fast and slow modes can be grouped and distinguished. In principle, we can 
express the right-hand member in ctj’s so that 

d N 

3^ = £/’«.• (7.3) 

» = 1 • 

Its success would lie in finding a way to decouple the modes through f 1 so that the amplitude of the 
uncoupled mode may evolve with its own characteristic time. The numerical procedure must deal 
with the convergence problem in its search for the normal mode of interest, apart from finding a 
robust method to identify the basis vector set which will not stay fixed and is highly temperature 
sensitive. The approach has been applied to methane combustion and also to dissociating hypersonic 
flows; it also has proven helpful in identifying the more important reactions in a complex, large 
system (Goussis et al. 1990, Gnoffo 1990). 


8. RAREFIED HYPERSONIC FLOW AND CONTINUUM EXTENSION 
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Owing to space limitation, the following will limit the discussion on rarefied flow to a minimum, and 
give preference to the issues on the continumm extension pertaining to fluid dynamic interest. 
Broader perspectives on RGD have been offered in earlier volumes of the Annual Review (Sherman 
1969, Kogan 1973, Bird 1978, Muntz 1989) and more recent works and reviews can be found in the 
1989 RDG Proceedings (Beylich 1990) and the forthcoming 1992 RGD Proceedings. 

8.1 DSMC as a Predictive Tool 

preliminary remarks The basic ideas of the DSMC method and numerical procedures were 
described in Bird’s (1976) monograph, and the underlying idea elucidated in several extensive reviews 
(Bird 1978, 1979, 1985, 1989). While the statistical errors present in a solution are expected to be 
proportional to the square root of the total number of simulated particles, i.e., V7T, an essential feature 
of the DSMC procedure is that the computation work is proportional to only the first power of N. 
The computer resource for a 2-D analysis is generally manageable in many institutes and universities, 
but is still much larger than that required for a NS calculation. Variants of this method have been 
developed for improved computational performance with vector-and parallel-processing computer 
architectures (Furlani & Lordi 1989, Baganoff & McDonald 1990, Boyd 1991, Wilmoth 1991, 
McDonald 1991). 


variable-hard -sphere model and inelastic collisions A novel feature which greatly increases 
the computing speed of Bird’s DSMC program is the variable-hard-sphere (VHS) model which treats 
the interaction of the molecules as a collision between two rigid spheres. This allows the post- 
collision relative velocity to be sampled from a uniform distribution in solid angle while the sphere 
radius is allowed to vary with the relative velocity c r so as to preserve the correct viscosity- 
temperature relation in translational equilibrium, i.e., oocc r ' 2lt/ . 

Extensions to molecules with rotational and vibration excitations have been made with the 
Borgnakke-Larsen (BL) (1075) phenomenological model, in which a fraction of the collisions is 
assumed elastic and the remainder inelastic. New values of the internal and translational energies are 
then sampled from an equilibrium distribution to be taken as the post-collision particle properties 
pertaining to the inelastic fraction. The characteristic collision number, Z, for the internal energy 
relaxation is roughly taken to be the reciprocal of the inelastic-collision fraction and may be chosen 
to match the relaxation -time estimated from experiment. The collision number Z has been commonly 
taken to be 5 for rotation and 50 for vibration. In their DSMC calculation for a shock front, Olynick 
et al. (1991) took into account the temperature and pressure dependence on the relaxation collision 
number in accordance with the rates used in the continuum model. Additionally, Boyd (1990) derived 
a velocity-dependent collision number which gives results in accordance with those of Parker (1959). 

Bypassing the BL phenomenological model, Boyd (1991) proposed a vibrational relaxation model 
for VHS application by modifying the Landau-Teller theory. It succeeds in determining vibrational 
cross sections for (one-step) activation and deactivation, which turn out to be independent of the 
intermolecular collision model. The model compared well with one experiment and represents a more 
detailed predicting method than the BL phenomenological model. An alternative modification of the 
BL model was offered by Chung et al. (1991) for a mixture which assumes a combining rule for the 
effective collision cross section. 

chemical reaction and radiation A form of collision theory of chemical physics that is consistent 
with the VHS model was used to convert the temperature-dependent rate constants to collision-energy 
dependent reaction cross sections Ambiguity in modelling the three-body collision arises. The latter, 
however, may not be of great concern for DSMC applications since the binary- scaling law prevails 
in most cases of nonequilibrium flow of interest (Gibson & Marrone 1962, Hall et al. 1962, Anderson 
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1989). A BL-type equilibrium distribution for the reaction products was again assumed; whether the 
model may allow vibrational-temperature dependence of the reaction rate is unclear. 

A similar idea was adopted in an extension of the DSMC to nonequilibrium flows with radiation 
(Bird 1987, Moss et al. 1988). Apart from issues on the BL phenomenological model, the task of 
tracking the large number of electronic states and of adjusting the constants for each rate to match 
data from available experiments and other continuum-based sets would appear insurmountable. 
Nevertheless, this task was accomplished and documented in detail. As a validation. Bird compared 
the spectral distribution of emitted radiation determined from an 1 1 -species air-chemistry calculation 
for an example corresponding to an AVCO experiment shown earlier in Figs. 9a,b, in which the thick 
dashes represent Bird’s (1987) results and appear to capture the trends well even in the nonequilibrium 
flow region. Inasmuch as the number of adjustable constants in the program is huge, an extensive 
sensitivity study is in order. Such an investigation was performed in part in Carlson & Hassan’s 
(1991) study where a scheme was introduced to determine the relaxation collision number "Z e " for the 
electronic-state excitation, which may reduce the degree of empiricism in the existing DSMC 
radiation model. 

Bird’s DSMC radiation model was applied by Moss et al. (1988) to predict the history of radiative 
and convective heating on an AFE vehicle which was investigated earlier by Hamilton et al. (1991) 
using various versions of continuum models (cf. §6) and corresponds to altitude 78-90. km and velocity 
7. 6-9.9 km/sec. The study indicates that the radiative heating is negligible compared to the 
convective heating; it becomes noticeable but still small during the peak heating period, where the 
stagnation convective and radiative heat fluxes are 0.19-0.21 MW/m 2 and .03-.04 MW/m , 
respectively, which are lower than the corresponding continuum VSL predictions by a factor of 0.40. 

8.2 More DSMC Calculations and Comparisons With Experiments 

A number of comparisons of recent DSMC calculations with experiments were reported in the 
Proceedings of the 17th International Symposium of Rarefied Gas Dynamics (Beylich 1990). There 
are several notable examples of comparisons from more recent studies. 

rarefied flow in nozzles The steady expansion of nitrogen from a 20° mozzle to a near- vacuum 
(with throat Knudsen number 2.3xl0" s ) was investigated experimentally and computationally using 
NS and DSMC (Boyd et al. 1991). Consistently good comparisons were found between DSMC 
calculations and Pitot pressure and flow angle measurements in and out of the nozzle. Whereas the 
solution is sensitive to the surface-interaction model, the fully diffuse wall model appears to be quite 
satisfactory. 

shock interference analyzed by DSMC The severe interference heating problem at a scramjet 
cowl lip noted in §3 (Fig. 2), the Edney Type IV problem, was attacked by Carlson & Wilmoth (1992) 
with DSMC calculations using 400,000 simulated molecules. The peak heating rate obtained appeared 
to be considerably lower than the experimental value q/q D ~ 30 (Wieting 1990). As the author noted, 
the grid/cell system employed may not be adequate, and a full NS calculation should have been made 
in this case. 

sphere dbag and wake structure DSMC calculations were made by Dogra et al. (1992) to 
compare with hypersonic sphere drags at low Knudsen numbers (M m = 11-13, Kn * 0.01-0.09) 
measured by Legg & Koppen waller (1970). The adequacy of the cell size and DSMC’s capability for 
describing unsteady separation are open to question. 

plume interaction, delta wing, afe vehiclb DSMC calculations and experiments on plume- 
freestream interaction were made by Campbell (1991); comparison in density distribution showed 
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qualitative agreement. A 3-D version of DSMC calculations for a delta wing (Celenligil& Moss 1991) 
were compared to a DLR wind-tunnel experiment at Mach 8.9, Knudsen number 0.02-2.0. Good 
agreement was found in lift, drag and surface heating rate; the surface temperature is believed to be 
near the stagnation value in this case. The aerodynamics of a "viscous optimized waverider" (cf. §2) 
in the rarefied gas dynamic regime was examined with DSMC calculations by Rault (1992). Using 
the F3 program of Bird (1990a,b), Rault gave a L/D of 0.24 at Mach 25, 100 km altitude (Kn = 0.01). 
It appears to be inferior to the delta wing last mentioned which has a L/D better than 0.5 at Kn less 
than 0.10. Interestingly, a 3-D DSMC calculation made recently for the very blunt AFE configuration 
(Celenligil et al. 1991) gave L/D = 0.212 at 100 km, not far from Rault’s value 0.24. 

8.3 More Detailed Validation of DSMC 

A more critical assessment of the DSMC method is to sample for the velocity-distribution functions 
f(u,v,w) and compare them with corresponding experimental data. Several sets of unpublished data 
for argon and helium ideal for this purpose were obtained earlier by E.P. Muntz for the partially 
integrated / within the shock-transition zone (inferred from measured intensity profiles of 
predominantly Doppler-broadened emission lines excited by an electron beam, known as electron 
fluorescence technique). Difficulties were encountered, however, in identifying the precise location 
at which each measurement was made, owing to a number of uncertainties related to the instrument 
and to the flow field nonuniformity. Thus the comparison study (Erwin et al. 1991; Pham-Van-Diep 
et al. 1989, 1991) served as a validation of the experiment procedure as well, especially since a 
convolving calculation procedure was adopted in some cases to identify the location that b est fits a 
particular experimentally determined (parallel or perpendicular) distribution function. In their 
analysis, Erwin et al. used differential cross sections based on Maitland -Smith (1981) potential, 
replacing the VHS collision model, which fits experimental viscosity data slightly better than other 
forms do. 

Two sets of predicted and experimental velocity distribution functions for a Mach 25 shock in 
helium are reproduced from Pham-Van-Diep et al. (1989) and shown in Figs. 13a,b, where the 
convolved parallel and perpendicular distributions 

/||(«) = J J fdvdw, / l («) = J j fdudw (8.1) 

are drawn in solid curves and dashes, respectively, with the corresponding experimental data shown 
as open circles and open triangles. The data set of Fig. 13a were identified with a location where the 
number density ratio a (n-n 1 )/n 2 -n 1 ) is 0.285, and the set of Fig. 13b pertaining to a further 
downstream station where n * 0.565. Similarly detailed agreement was found with helium at Mach 
1.59 and argon at Mach 7.18 (Erwin et al. 1991). These close comparisons indicate the remarkable 
ability of DSMC to predict population of scattered atoms in this highly nonequilibrium state; they also 
reveal unmistakably the Mott-Smith (1951) type bimodal-like distribution, which signifies the 
persistent influence of the upstream and downstream states. It remains to be seen if the VHS version 
of the DSMC may also produce similarly encouraging comparison. [A preliminary study with VHS 
(Muntz et al. 1991) indicated general agreement with noticeable differences in the vicinity of the fi 
maximum.] 

8.4 Continuum Extension to Rarefied Hypersonic Flow 

The foregoing examples have demonstrated that the DSMC method can treat problems normally 
handled by the NS- based equations, but it demands large and costly computer resources. For 
example, one of Carlson et al. (1992) DSMC calculations took 33 days on a dedicated Sun SPARC 
station-2, and Celenligil's et al. (1991b) 3-D calculation needed 35 CPU hours on CRAY-2. The NS 
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Fig. 13 Predicted and experimental velocity distribu- 
tions (from Pham-Van-Diep et. al. 1989): (a) 
h = 0.285, (b) n = 0.565; experiment (symbols), 
DSMC parallel (solid curve), DSMC perpendicular 
(dashes). 
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Fig. 14 Example comparing DSMC, VSL, FVSL, and 
NS calculations: (a) overal temperature distribu- 
tions along stagnation streamline with Rs = 13m, 
Uoo — 7.5km/sec., alt. 92km, T w ~ lOOO^i: , non- 
catalytic surface; DSMC (open circle), VSL (solid 
curve), FVSL (dashes), [reproduced from Moss L 
Bird (1985)] (b) Surface heat-transfer coefficients 
as function of Knudsen number for f/oo = 7.2 - 
7.5km/sec., alt. 75 - 115km; DSMC (open circle), 
VSL without slip (squares), NS without slip (trian- 
gles with slashes), NS with slip (triangles) (repro- 
duced from Gupta & Sunmonds 1986). 
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based solutions/approximations have proven usef ul in low-density hypersonic flow studies (e.g. Cheng 
1966) and compared reasonably well with surface heat flux measurements in hypersonic flows 
(Wittiliff & Wilson 1962, Vidal et al. 1963). With DSMC calculations, it becomes possible to assess 
the NS based results and other continuum extensions and identify their applicability domains. This 
was, in fact, accomplished by Moss & Bird (1985), Gupta & Simmonds (1986) Moss et al. (1987) and 
many subsequent workers, comparing DSMC, NS and NS- based viscous shock layer (VSL) calculations 
for the blunt nose region. Certain notions and concepts used in the subsequent discussions need to 
be clarified. 

viscous and fully- viscous shock layess Just as the inviscid shock layer in the classical theory, 
the concept of a thin shock layer is applicable to the viscous, heat-conducting flow region between 
the shock and the surface of a blunt/non-slender body, provided the density level there is much 
higher than that upstream. Viscous formulations based on thin-shock layer approximations could all 
be called viscous shock layer, but a distinction must be made between the version which use the 
viscous modified Rankine-Hugoniot relation at the outer boundary and that which assumes the 
inviscid shock relation. The latter version has been called the viscous shock layer (VSL) by Moss and 
coworkers, although being somewhat inconsistent in a strict sense. Thus the former version with the 
viscous modified shock condition will be referred to as the fully viscous shock layer (FVSL). 

the shock slips Owing to its resemblance to the wall slip, the change in the shock boundary 
conditions in the FVSL formulation, which includes corrections in tangential velocity and total 
enthalpy, has been called by Davis (1970), Moss and others "shock slips". Underlying these modified 
shock conditions is the stipulation that the density in the shock interior (structure) is low and 
comparable to the free-stream level; thus the tangential components of the mass, momentum and 
energy fluxes can little affect the balances in the normal flux components, as long as the thickness 
of the shock (structure) is small compared to the shock/body radius of curvature (even if the shock 
thickness becomes comparable to the thickness of the shock layer (Cheng 1961)]. This also implies 
that a shock -capturing NS solution should provide the shock slips correctly, even if the shock 
structure so obtained may not be physically correct on a kinetic -theory basis. The FVSL with the 
shock slips will provide a framework where a kinetic-theory base for the continuum extension can 
be found. 

strong wall cooling and wall slips A surface with low wall-to-stagnation temperature ratio 
(T w /T 0 «l) is more of relevance to hypersonic flight than one with a nearly insulated wall; the strong 
cooling also makes the thin shock-layer analyses a better approximation for the FVSL, since the 
shock-layer density level is raised substantially through wall cooling. An additional consequence of 
strong wall cooling is in the limiting of the wall-slip influence on the FVSL and VSL, to a relative 
order of (Cheng 1966) 


where e in above is the small parameter in theWock-Iayer theory denoting the ratio of the free stream 
density to a typical shock-layer density. 

comparison of ns-based and DSMC calculations To indicate the degree to, and the manner in 
which the NS-based predictions may differ from one another and from DSMC calculations in regime 
where departure from local translational equilibrium is significant. Figs. 14a,b, reproduce from Moss 
& Bird (1985) a comparison for temperature distribution along a stagnation streamline, and from 
Gupta & Simmond (1986) a comparison for the stagnation-point heat-transfer coefficient. The 
example of Fig. 14a has a strongly cooled noncatalytic nose with a 1.3 meter radius, at 7.5 km/sec 
speed, 92.4 km altitude. Here, Moss’ NS-based VSL predictions (in solid lines) gives noticeably 
higher overall (effective) temperature than the DSMC data (in open circles) at the outer part of the 
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VLS; to indicate the shock-slip effect, a corresponding FVSL profile based on a binary reacting gas 
model (Cheng 1963) is also included as dashes; the viscous and heat-conduction corrections at the 
outer edge decrease the temperature in the outer portion of the FVSL and reduce its thickness. This 
figure serves to illustrate the great extent of the shock structure in this type of flow, which may be 
considered being almost twice the shock-layer thickness in this case. But owing to its lower density 
level, the influence on the downstream flow is limited. Figure 14b compares the heat-transfer 
coefficient C H over a wide range of Knudsen numbers for a nose radius about 1.36 meter at speed 
7.5 km/sec. Comparing with the DSMC results (in open circles), the VSL without wall slips and shock 
slips (in open squares) begins to depart at a Kn - 0.025, as Kn increases. The NS without wall slips 
(in open triangles) remains close to the main trend up to Kn ~ 0.25 where the shock slips are expected 
to be significant; this is consistent with the minor wall-slip effect (8.1) observed earlier. Surprisingly, 
the solution for NS with wall slips agrees with DSMC all the way up to Kn - 1. However, in a similar 
study by Lee et al. (1990), the NS with wall slips reaches the limit C H * 1 at Kn » 2.5 and suggests 
a peculiar trend of overshooting the limit at Kn > 2.5. 


8.5 On Kinetic-Theory Basis of NS, Burnett and Thirteen-Moment Equations 

Early research on improvement of the NS description for rarefied flows by applying the Burnett 
(1936) and Grad’s (1949) thirteen moment equations was not successful, and the general belief was 
that these higher-order equations from kinetic theory do not seem able to predict when the NS 
relations break down (see Schaaf & Chambr£ 1958, p. 718). This perception has apparently changed 
by recent studies with the help of modern CFD and with DSMC calculations, notably the study of the 
plane shock structure by Fiscko & Chapman (1988a,b). 

series expansion near continuum limit The gas-kinetic base of the NS and Burnett equations is 
the formalism of the Chapman-Enskog expansion and the higher-order development of the velocity- 
distribution function for a short particle-collision time, A/c or /i/p, compared to the flow 
characteristic time (Chapman & Cowling 1953, Vincenti & Kruger 1965, Kogan 1969, Ferziger & 
Kaper 1972). 


thirteen-moment system Grad’s system of 13-rooment equations is a particular set of moments 
of the Boltzmann equations— the Maxwell transfer equations, in which closure is achieved with the 
help of a form of the distribution function 

/ = p e - ^ c ’[l + A ijCiCj + B iCi + C<c 2 c,-] (8.2a) 

where ft - (2RT)' 1 , C; is a thermal velocity component, c 2 = c^+Cj^Cj 2 , and the polynomial 
coefficients can be identified with stress-tensor components and heat fluxes as 


» _ Pij n ,_ £<_ c =-£ — 

iJ ~2 pRT' pRT' pRTSRT' 


(8.2b) 


The polynomial inside the square bracket came as a truncated Hermite polynomial expansion carried 
out for the twentieth moment in Grad’s (1949) original work.' Equation (8.2) is precisely the form for 
f needed in the Chapman-Enskog theory for the derivation of the NS and Fourier constitutive 
relations and the evaluation of the viscosity and heat-conductivity coefficients. Whereas, the Pjj, q,, 
Uj and p are to be solvewd as unknowns in the PDE’s of the 13-moment system, which does not require 
/i/p or Kn - A/L to be small as long as (8.2) remains adequate. This amounts to allowing nonvanishing 


Pij/P, qi/'/iRT - 0(1) 


(8.3) 


Interestingly, the full Burnett equations, including the boundary conditions for velocity and 



temperature, are contained in (derivable from) Grad’s 13-moment equations in the case of a Maxwell 
gas for asymptotically small p/p or Ip^/pl, even though Burnett’s original theory involved a form 
far more complicated and of order higher than that in (8.2), as noted by Schaaf & Chambre (1958) 
and made more explicit by Yang (1992). The corrections for the non-Maxwell monatomic gas in both 
Burnett and 13-moment equations are rather small (Chapman & Cowling 1953, Fiscko & Chapman 
1988; Chen, private communication). 

BASIC issues on burnett and 13-moment EQUATIONS The Burnett equations as well as the extension 
to the Super-Burnett equations (Simon 1976) are valid only as successive corrections to the Euler 
equations for a nearly inviscid flow (outside of boundary and shear layers), inasmuch as 
inviscid/isentropic relations were used to simplified terms Dp^/Dt and Dq/Dt in the final form of 
the constitutive relations; it is not strictly applicable to either a fully viscous or boundary- layer 
region. This makes it more restrictive than merely requiring Ip^/pl being small. More critical is the 
unresolved issue with the proper boundary conditions when the Burnett system is solved as full 
equations without further approximations for Kn * 0. As Schaaf & Chambre (1958) noted, at the 
Burnett level, terms of an order higher in derivatie occur in each equation and an additional boundary 
condition must be prescribed. Thus a nonuniqueness problem will arise unless the system possesses 
a special property for the exception. The issue cannot be settled by simply demonstrating the 
existence of a solution. 

For application to supersonic flows. Grad’s 13-moment equations fail to yield a normal shock 
structure at Mach number exceeding 1.65, as is well known (Grad 1952). The system without 
modification cannot provide a base for the analysis of the entire flow field which include the shock 
structure. The discussion in §8.7 will examine its applicability to the study of flow behind shock. 

8.6 Burnett Equations as CFD Model for Rarefied Hypersonic Flows 

NS calculations have been known to give poor shock-structure descriptions when compared to 
particle-simulation and experimental results (Schaaf & Chambre 1958, Bird 1978, Muntz 1989). 
Using flux-splitting technique, Fiscko & Chapman (1988a,b) found that Burnett equations provide 
much greater accuracy than NS equations for 1-D shock structure in a monatomic gas, although its 
degree of improvement over NS varies, depending on the Mach number and the viscosity-temperature 
law of the gas, and on the flow quantities of interest. Figure 15 reproduces from Fiscko & Chapman 
(1988a) a comparison of the inverse density thicknesses (IDT) for an argon shock predicted by Burnett 
equations (in solid curves), NS equations (in fine dots) and DSMC (in open circles) over a wide range 
of Mach numbers (cf. Muntz 1989 for definition of IDT). The encouraging agreement with the 
particle-simulation results suggests that a CFD model based on the Burnett equations may provide a 
much improved shock-capturing capability in a rarefied hypersonic flow where a realistic description 
of the shock structure at high temperature is of vital aerothermodynamic interest. Apart from the 
issue on the surface boundary conditions, there is another obstacle to this extension, noted and 
resolved by Zhong et al. (1991). Namely, the numerical solution is linearly unstable to disturbances 
with wave length comparable, or less than, the mean free path (Bobylev 1982, Foch 1973), and this 
difficulty was believed to place a handicap on Fiscko’s eaTlier calculations; it would also pose a 
computation problem if one desires to study the higher altitude effects for the same body geometry 
using a fixed grid. Zhong et al. overcame this problem by adding several higher-order derivative 
terms to the Burnett constitutive equations. These added terms have forms similar to certain terms 
at the Super-Burnett level but with different coefficients (and signs). The augmented Burnett 
constitutive equations/relations for the deviatorial stress and heat-flux components are written as 
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Fig. 16 Comparison of stabilized Burnett, Navier- 
Stokes and particle simulation calculations along 
stagnation streamline for case IV — 25, a = 
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where supercripts (1) and (2) refer to the NS and Burnett levels, respectively, and (a) indicates the 
added terms. In the 1-D case, for example, o u W and q^ 2 ) contain terms like 
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and the added terms are 
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pp L p J 


The set w 7 = 2/9, 0 6 =-5/8 and 0 7 = 11/16, which is by no means the unique choice, proves to make 
the linearized Burnett equations stable. Incidentally, this set turns out to be precisely that given by 
Wang-Chang (1948) for a Maxwell gas, of which the 0 7 value did not agree with the correct value 
-157/1 16. The augmented Burnett results enjoy a much higher accuracy than Fiscko’s earlier version 
and are marked also in Fig. 15 as open triangles. 


More interesting are perhaps the 2-D examples in Zhong et al. (1991) which compares the NS, 
augmented Burnett, and particle-simulation calculations. They are among the first Burnett solutions 
to hypersonic blunt-body problems to appear, and to which the boundary-condition issues must be 
addressed. The (M^Kn) pairs considered in the four cases are (4,067x1 0" 4 )( 1 0,0. 1 0)( 10,1 .2) and 
(25,0.28), referred to as cases I, II, III and IV, respectively. A constant specific heat ratio 7 = 1 .40 was 
assumed, except in case IV where rotational relaxation is allowed in both Burnett and the particle 
simulations; whether the relaxation models can be identified to those in Lumpkin & Chapman (1991 ) 
and Lumpkin et al. (1989) is unclear. For these examples, the surface temperature is not low, being 
in the range T w /T c ~ 1/3 - 1/2. Here, the Burnett and NS solutions shared the same slip boundary 
conditions, while the corresponding surface conditions for the particle simulation calculation was not 
explained. Of interest are the density and temperature distributions along the stagnation streamline 
in cases II, III and IV where the differences among the three solution sets are surprisingly small, 
especially in the density profiles. Even in the thick shock interior, where discrepancies of the NS 
temperature predictions are noticeable, the NS solution actually describes the temperature profile not 
so poorly. 

Figure 16 reproduced from Zhong et al. (1991) compares the density and temperature for case IV. 
Of vital engineering interest is that temperature profiles of the three solutions sets become 
indistinguishable from one another downstream of the shock structure in all cases. This is rather 
surprising, not only because NS results have never before been shown being in agreement so well with 
Burnett and particle simulation in the shock interior, but also because the FVSL parameter K . 2 (Cheng 
1966, see below) is of unit order or smaller in cases II-IV, which signifies large departures from 
translational equilibrium. These could have resulted from the Sutherland viscosity law assumed or 
from the particular iterative procedure used which avoid the additional boundary conditions by 
extrpolating from the interior; it remains to be further examined. While the study shows Burnett 
equations can yield a solution closely matched the particle-simulation results, the same study also 
indicates that NS calculations can predict the shock and flow structures almost as well as the Burnetts. 
This interesting CFD work falls short in settling the obvious uniqueness issue which cannot be 
answered by simply demonstrating the solution’s existence, as noted earlier. For the steady Couette 
flow problem, Lee (private communication) shows that, even with the extrapolation technique, the 
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solution is not unique and depends on the initial input for the iteration. 


8.7 Thirteen-Moment Equations as Basis for Viscous Shock Layer Theory 


need or gas kinetic basis As noted earlier, the use of FVSL, VSL as well as the full NS equations 
in rarefied hypersonic flow analyses anticipates the viscous and other molecular transport effects to 
rank equally with the convective processes, which also implies condition (8.3). The very same 
condition indicates, however, a large departure from the translational equilibrium, invalidating the 
gas-kinetic base for the NS and the Burnett equations. Grad’s 13-moment theory, which allows 
condition (8.3), may therefore serve a kinetic-theory base for viscous shock-layer anslyses. Another 
reason for using the 13-moment equations is the absence of ambiguity in the type and number of 
admissible boundary conditions on a body surface. By considering the number of characteristics 
reaching the boundary. Grad (1949) showed that the number and type in question are the same as in 
NS; this can be confirmed by examining the nature of the 13-moment equations in a Couette flow 
[Cheng et al. 1989, Eq. (4.35)]. Using scales typical of the flow in a shock layer, it is possible to 
express the order of magnitude of Ip^/pl more explicitly in the form of the reciprocal of a local 
Reynolds number 


Poo Uoo X 
X = f 

Po 



(8.5) 


where 8* is the shock or body incidence angle, the subscripts refer to a suitable reference condition 
and x is a distance or reference length. With x replaced by the nose radius R N , X is identified with 
the K 2 in Cheng’s (1961,1963,1964) early NS-based theory. Thus X or K 2 directly control both the 
shock slip and departure from translational equilibrium. The Knudsen number Kn s A^/L an ^ 

V s x/M^ commonly used in rarefied gas dynamics are related to T as 

(x) -1 , V ~ (xco®^) -1 / 2 ^ 8 ' 6 ^ 


the thin-layer approximation and QUASi-i-D shock STRUCTURE Using scales appropriate to the 
shock structure, thin-layer approximations can be applied to derive the equations governing the 
quasi- 1-D shock structure. They can be integrated to arrive at the modified Rankine-Hugoniot 
relation allowing the shear-stress and heat-flux contributions immediately behind the shock, 
identified as the "shock slips". Using conventional notations in the shock-layer theory, with y and 
v referring to the coordinate and velocity component in the direction normal to the body surface. The 
latter conditions read 


u - Uoo = pia/mi, w-wi=p3t/m l (8.7) 

H — Hoc = (upu + wps2 + q2)/ m \ 


where mj* p 00 v 1 is the component of free-stream velocity normal to the surface, q 2 is the normal heat 
flux component, p 12 and p 32 are pressure tensor components associated with (x,y) and (z,y) 
respectively. Subject to error of order e, this provides the outer boundary conditions for the shock- 
layer flow, irrespective of the gas-kinetic model in the shock interior (Cheng et al. 1991). The 
remaining shock condition is P 22 * nijVj, where the normal pressure-tensor component P 22 is not the 
thermodynamic pressure p, owing to translational noequilibrium. 

fvsl based on 1S-MOMENT equations The governing equations for the shock layer are, to the 
leading order, essentially the same as in viscous shock-layer theories, except that constitutive relations 
expressing p 12 , p 32 , p 22 and q 2 in terms of flow gradients must now be replaced by the 13-moment 
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leading order, essentially the same as in viscous shock-layer theories, except that constitutive relations 
expressing p 12 , p 32 , p 22 and q 2 in terms of flow gradients must now be replaced by the 13-moment 
theory. Under the thin-layer approximation and with the formalism of a small e familiar from the 
shock-layer theory, i.e., considering R/C p being small, the constitutive relations in question can be 
reduced simply to (Cheng et al. 1989, 1991) 


Pi 2 = 


Pn du 
P dy 


P32 = — 


P22 dw 
P ^ dy 


9 2 = 


P 22 ,dT 
P dy 


(8.8a,b) 

(8.8c) 
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The above expressions differ from the corresponding NS based relations mainly in the appearance of 
the common factor P 22 /p to be determined through a nonlinear relation to the velocity and 
temperature gradients (8.8d). It is, in fact, through the p, not P 22 , that translational nonequilibrium 
will affect the dynamics and thermodynamics of a viscous shock layer. Except for the last equation, 
the above constitutive relations involve derivatives of u and T no higher than the first order and 
differ from those in the Burnett theory. We note that the last (8.8d) is a thermal stress (Kogan 1969) 
and is formally a higher order term, but is included for its exceptionally large coefficient and physical 
significance. To complete the formulation, we may apply wall-slip equations from Grad’s wall model 
suitably simplified in a manner consistent with (8.7)-(8.8). These wall-slip effects, however, can 
influence the shock-layer flow at the most of relative order {eT w /T^ 1 ^ as indicated earlier. They 
have been demonstrated to be negligibly small even for exceptionally large accommodation 
coefficients (Cheng & Wong 1988). 

The observation made above on the nonequilibrium influence through the ratio P 22 /p indicates the 
possibility of successfully correlating a kinetic-based shock-layer flow with a NS-based flow. Its key 
lies in the recognition that the reciprocal of the density p or p always appears in a product together 
with normal derivative d/dy both in the governing PDE’s and in the outer boundary conditions 
involving the shock slips. The e or p may then be eliminated in this case by the use of Dorodnitsyn 
or von Mises variables, through which the governing equation system is transformed to a NS-based 
system, with slightly different wall-slip boundary conditions. The latter is inconsequential for the 
strongly cooled surface of interest. This correlation principle holds for the tangential velocity 
components and the enthalpy, and also for the major stress components and the heat flux. A 
consequence is that the skin friction and surface heating rate are predictable from the NS-based 
equations, being unaffected by the translational nonequilibrium to the leading order, although the 
streamline pattern and shock-layer thickness will be accordingly displaced. The present version of 
the theory therefore provides a kinetic-theory base for explaining the good agreement of the NS- 
based FVSL analyses with early heat-transfer measurements. 

FLAT PLATE AT INCIDENCE: CORRELATING DSMC AND NS-BASED CALCULATIONS The flow about a flat 
plate at 40° attack angle was studied by Cheng et al. (1989,1990,1991), and also Cheng (1989), as a 
generic lifting surface problem in rarefied hypersonic flow, for which FVSL, parabolized (thin-layer) 
NS, and time-accurate NS calculations were made and compared to corresponding DSMC 
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computations (Dogra & Moss 1989, Dogra et al. 1989) and examples include (monatomic) Maxwell gas 
and diatomic-gas model with 7 = 1.40 and Pr * 0.72 (which implies a fast rotation- translation energy 
transfer), and cases with different viscosity- temperature relations, wall temperatures, and wall-slip 
models. While some of the DSMC data were generated from a code which allows internal excitation 
and chemical reaction of a model air, the latter effects on aerothermodynamics are negligible even 
at flow speed of 7.5 km/sec, owing to the high degree of gas rarefaction at altitudes 90-130 km. The 
coefficients of heat transfer and skin friction, C H and C { , from three sets of DSMC data and six sets 
of NS-based results, each with different M^, Re and viscosity law were determined over a wide range 
of X, 0.2- 1 0 (X = 1 0" s - 1 0, if the DSMC data near the collision-free limit is also included). Figure 1 7 
presents the correlation/comparison in C H as function of X; even though a perfect Reynolds analogy 
is not expected, the corresponding correlation for C f in this case turns out to be almost 
indistinguishable from that in Fig. 17 (except a portion of DSMC calculation for the 100 km altitude 
which defies explanation) 

The lift-to-drag ratio of the plate at 40° incidence computed from the integrated normal and 
tangential forces on the windward side are reproduced from Cheng et al. (1991) as a function of X in 
solid curves, which agrees exceedingly well with Dogra et al. (1989) DSMC calculation for a one- 
meter plate (in open circles with slashes) over the entire X-range (2xl0' 2 through 15). Also included 
as a solid curve is the L/D computed for a plate at 20° attack angle. Of interest are L/D values 
computed for ten widely different planforms based on the 2-D distributions of skin friction and 
normal force (reproduced from Cheng (1991). These calculations (Hoover et al. 1992) used a version 
of the strip method which is the result of the 3-D FV SL theory for a flat-bottom surface (Cheng et 
al. 1991). When the distance x in the 2-D problem were taken to be the span-averaged chord (SAC), 
the ten L/D values fall into the vicinity of the strictly 2-D results. This is surprising, because as 
altitude or Knudsen number increases, the planform dimension and shape are expected to strongly 
affect the skin friction hence the L/D. This insensitivity of the L/D on planform indicates a way to 
identify the bridging function for planar lifting surface (Warr 1970, Wilhite et al. 1985, Potter 1988), 
which is provided here by the 2-D data for an inclined plate. 

The flows over an aligned flat plate and on the lee side of a flat surface remain as examples for 
which still more can be learned from critical examination of the continuun-extension and particle 
simulation calculations. 


9. CONCLUDING REMARKS 


In this article, an attempt has been made to reflect on current focuses in certain areas of hypersonic 
flow research by examining the recent works and their issues. The study addresses aspects of viscous 
interaction, flow instability, and nonequilibrium aerothermodynamics pertaining to theoretical 
interest. The field is a diverse one, and many exciting works may have either escaped the writer’s 
notice or been abandoned for the sake of space. As noted in the text, students of hypersonic viscous 
flow must face the transition problems towards the two opposite ends of the Reynolds or Knudsen 
number range (§§5,8), which represents two regimes where unresolved fluid/gas dynamic problems 
abound. Central to the hypersonic flow studies is high -temperature physical gas dynamics (§6); here, 
a number of issues on modelling the intermolecular potentials and inelastic collisions remain the 
obstacles to quantitative predictions. Research in combustion and scram jet propulsion will certainly 
be benefitted by advances in turbulent mixing and new CFD strategies on multi-scaled complex 
reactions. 

Apart from many omissions, the view and interest expressed by the author has been limited to the 
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theorectical side. However, even for the sake of theoretical development, the lack of pertinent 
experimental data in the right energy and density ranges is believed to be among the major obstacles 
to progress in aerothermodynamic research for hypersonic flight. Flight experiments performed by 
Space-Shuttle missions and the AFE Program referred to often in the text have been extremely 
valuable. To enable laboratory simulation of nonequilibrium effects anticipated for transatmospheric 
flight, facilities capable to generate high enthalpy flow at density levels higher than in existing 
laboratories are needed (Hornung 1988). A new free-piston shock tunnel capable to realize a test- 
section stagnation temperature of 10 6 K° at Reynolds number 50xl0 6 /cm art ‘being completed and 
preliminary tests hatfbegun (H. Hornung et al. 1992). Another laboratory study worthy of note as well 
as theoretical support is the nonequilibrium flow experiment of iodine vapor which has low activation 
energies for vibrational excitation and dissociation, and can be studied in a laboratory with modest 
resources (Pham-Van-Diep et al. 1992). , 
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